diff --git a/src/analysis/fba.jl b/src/analysis/fba.jl
index f86e9aec3c859dc9528ccbb3cb7609554726bec0..ad8a5c74887958d965d0f67f3678f6007be06bc1 100644
--- a/src/analysis/fba.jl
+++ b/src/analysis/fba.jl
@@ -22,10 +22,9 @@ Arguments are passed to [`flux_balance_analysis`](@ref).
 """
 function flux_balance_analysis_vec(args...)::Union{Vector{Float64},Nothing}
     optmodel = flux_balance_analysis(args...)
-    vars = optmodel[:x]
-
+    
     termination_status(optmodel) in [MOI.OPTIMAL, MOI.LOCALLY_SOLVED] || return nothing
-    value.(vars)
+    value.(optmodel[:x])
 end
 
 """
@@ -69,7 +68,7 @@ function flux_balance_analysis(
 )
     # get core optimization problem
     cbm = make_optimization_model(model, optimizer)
-
+    
     # apply callbacks - user can also just put in a function
     if typeof(modifications) <: Vector
         for mod in modifications
@@ -80,16 +79,22 @@ function flux_balance_analysis(
     end
 
     optimize!(cbm)
+    
+    return cbm
+end
 
-    status = (
-        termination_status(cbm) == MOI.OPTIMAL ||
-        termination_status(cbm) == MOI.LOCALLY_SOLVED
-    )
+function flux_balance_analysis_vec()
 
-    if status
-        return map_fluxes(v, model)
-    else
-        @warn "Optimization issues occurred."
-        return Dict{String,Float64}()
-    end
+end
+
+status = (
+    termination_status(cbm) == MOI.OPTIMAL ||
+    termination_status(cbm) == MOI.LOCALLY_SOLVED
+)
+
+if status
+    return map_fluxes(cbm[:v], model)
+else
+    @warn "Optimization issues occurred."
+    return Dict{String,Float64}()
 end
diff --git a/src/analysis/fva.jl b/src/analysis/fva.jl
index 380fac2880af3ae0759bfd5b078415cf90d675ce..9c739b3054292f47847ab096393181863782deb2 100644
--- a/src/analysis/fva.jl
+++ b/src/analysis/fva.jl
@@ -134,55 +134,20 @@ fva_max, fva_min = fva(model, biomass, optimizer; solver_attributes=atts)
 function fva(
     model::CobraModel,
     optimizer;
-    objective_func::Union{Reaction,Array{Reaction,1}} = Reaction[],
     optimum_bound = 0.9999,
-    weights = Float64[],
-    solver_attributes = Dict{Any,Any}(),
-    constraints = Dict{String,Tuple{Float64,Float64}}(),
-    sense = MOI.MAX_SENSE,
+    modifications = [(model, opt_model) -> nothing]
 )
-    cbm = make_optimization_model(model, optimizer, sense = sense)
+    # get core optimization problem
+    cbm = make_optimization_model(model, optimizer)
     v = cbm[:x]
 
-    if !isempty(solver_attributes) # set other attributes
-        for (k, v) in solver_attributes
-            set_optimizer_attribute(cbm, k, v)
-        end
-    end
-
-    # set additional constraints
-    for (rxnid, con) in constraints
-        ind = model.reactions[findfirst(model.reactions, rxnid)]
-        set_bound(ind, cbm; lb = con[1], ub = con[2])
-    end
-
-    # if an objective function is supplied, modify the default objective
-    if typeof(objective_func) == Reaction || !isempty(objective_func)
-        # ensure that an array of objective indices are fed in
-        if typeof(objective_func) == Reaction
-            objective_indices = [model[objective_func]]
-        else
-            objective_indices = [model[rxn] for rxn in objective_func]
-        end
-
-        if isempty(weights)
-            weights = ones(length(objective_indices))
-        end
-        opt_weights = zeros(length(model.reactions))
-
-        # update the objective function tracker
-        # don't update model objective function - silly thing to do
-        wcounter = 1
-        for i in eachindex(model.reactions)
-            if i in objective_indices
-                # model.reactions[i].objective_coefficient = weights[wcounter]
-                opt_weights[i] = weights[wcounter]
-                wcounter += 1
-                # else
-                # model.reactions[i].objective_coefficient = 0.0
-            end
-        end
-        @objective(cbm, sense, sum(opt_weights[i] * v[i] for i in objective_indices))
+    # apply callbacks - user can also just put in a function
+    if typeof(modifications) <: Vector
+        for mod in modifications
+            mod(model, cbm)
+        end    
+    else
+        modifications(model, cbm)
     end
 
     optimize!(cbm)