diff --git a/core/setParam.m b/core/setParam.m index e24ad2ca..734853b7 100755 --- a/core/setParam.m +++ b/core/setParam.m @@ -39,6 +39,9 @@ else rxnList=convertCharArray(rxnList); end +if isempty(rxnList) + return; +end %Allow to set several parameters to the same value if numel(rxnList)~=numel(params) && numel(params)~=1 @@ -59,23 +62,16 @@ %as we do not want to throw errors if matches fail indexes=zeros(numel(rxnList),1); -for i=1:numel(rxnList) - index=find(strcmp(rxnList{i},model.rxns),1); - if ~isempty(index) - indexes(i)=index; - else - indexes(i)=-1; - EM=['Reaction ' rxnList{i} ' is not present in the reaction list']; - dispEM(EM,false); - end +[Lia,Locb] = ismember(rxnList,model.rxns); +indexes(Lia) = Locb; +if any(~Lia) + params(~Lia)=[]; + indexes(~Lia)=[]; + paramType(~Lia)=[]; + dispEM('Reactions not present in model, will be ignored:',false,rxnLise(~Lia)); end -%Remove the reactions that were not found -params(indexes==-1)=[]; -indexes(indexes==-1)=[]; -paramType(indexes==-1)=[]; %Change the parameters - if ~isempty(indexes) if contains(paramType,'obj') model.c=zeros(numel(model.c),1); % parameter is changed, not added diff --git a/doc/core/setParam.html b/doc/core/setParam.html index 73b63b7e..96f0a107 100644 --- a/doc/core/setParam.html +++ b/doc/core/setParam.html @@ -111,96 +111,92 @@

SOURCE CODE ^else 0040 rxnList=convertCharArray(rxnList); 0041 end -0042 -0043 %Allow to set several parameters to the same value -0044 if numel(rxnList)~=numel(params) && numel(params)~=1 -0045 EM='The number of parameter values and the number of reactions must be the same'; -0046 dispEM(EM); -0047 end -0048 -0049 if length(rxnList)>1 -0050 if length(paramType)==1 -0051 paramType(1:length(rxnList))=paramType; -0052 end -0053 if length(params)==1 -0054 params(1:length(rxnList))=params; +0042 if isempty(rxnList) +0043 return; +0044 end +0045 +0046 %Allow to set several parameters to the same value +0047 if numel(rxnList)~=numel(params) && numel(params)~=1 +0048 EM='The number of parameter values and the number of reactions must be the same'; +0049 dispEM(EM); +0050 end +0051 +0052 if length(rxnList)>1 +0053 if length(paramType)==1 +0054 paramType(1:length(rxnList))=paramType; 0055 end -0056 end -0057 -0058 %Find the indexes for the reactions in rxnList. Do not use getIndexes -0059 %as we do not want to throw errors if matches fail -0060 indexes=zeros(numel(rxnList),1); -0061 -0062 for i=1:numel(rxnList) -0063 index=find(strcmp(rxnList{i},model.rxns),1); -0064 if ~isempty(index) -0065 indexes(i)=index; -0066 else -0067 indexes(i)=-1; -0068 EM=['Reaction ' rxnList{i} ' is not present in the reaction list']; -0069 dispEM(EM,false); -0070 end -0071 end -0072 -0073 %Remove the reactions that were not found -0074 params(indexes==-1)=[]; -0075 indexes(indexes==-1)=[]; -0076 paramType(indexes==-1)=[]; -0077 %Change the parameters -0078 -0079 if ~isempty(indexes) -0080 if contains(paramType,'obj') -0081 model.c=zeros(numel(model.c),1); % parameter is changed, not added -0082 end -0083 for j=1:length(paramType) -0084 if strcmpi(paramType{j},'eq') +0056 if length(params)==1 +0057 params(1:length(rxnList))=params; +0058 end +0059 end +0060 +0061 %Find the indexes for the reactions in rxnList. Do not use getIndexes +0062 %as we do not want to throw errors if matches fail +0063 indexes=zeros(numel(rxnList),1); +0064 +0065 [Lia,Locb] = ismember(rxnList,model.rxns); +0066 indexes(Lia) = Locb; +0067 if any(~Lia) +0068 params(~Lia)=[]; +0069 indexes(~Lia)=[]; +0070 paramType(~Lia)=[]; +0071 dispEM('Reactions not present in model, will be ignored:',false,rxnLise(~Lia)); +0072 end +0073 +0074 %Change the parameters +0075 if ~isempty(indexes) +0076 if contains(paramType,'obj') +0077 model.c=zeros(numel(model.c),1); % parameter is changed, not added +0078 end +0079 for j=1:length(paramType) +0080 if strcmpi(paramType{j},'eq') +0081 model.lb(indexes(j))=params(j); +0082 model.ub(indexes(j))=params(j); +0083 end +0084 if strcmpi(paramType{j},'lb') 0085 model.lb(indexes(j))=params(j); -0086 model.ub(indexes(j))=params(j); -0087 end -0088 if strcmpi(paramType{j},'lb') -0089 model.lb(indexes(j))=params(j); -0090 end -0091 if strcmpi(paramType{j},'ub') -0092 model.ub(indexes(j))=params(j); -0093 end -0094 if strcmpi(paramType{j},'obj') -0095 model.c(indexes(j))=params(j); +0086 end +0087 if strcmpi(paramType{j},'ub') +0088 model.ub(indexes(j))=params(j); +0089 end +0090 if strcmpi(paramType{j},'obj') +0091 model.c(indexes(j))=params(j); +0092 end +0093 if strcmpi(paramType{j},'rev') +0094 %Non-zero values are interpreted as reversible +0095 model.rev(indexes(j))=params(j)~=0; 0096 end -0097 if strcmpi(paramType{j},'rev') -0098 %Non-zero values are interpreted as reversible -0099 model.rev(indexes(j))=params(j)~=0; -0100 end -0101 if strcmpi(paramType{j},'var') -0102 if params(j) < 0 -0103 model.lb(indexes(j)) = params(j) * (1+var/200); -0104 model.ub(indexes(j)) = params(j) * (1-var/200); -0105 else -0106 model.lb(indexes(j)) = params(j) * (1-var/200); -0107 model.ub(indexes(j)) = params(j) * (1+var/200); -0108 end -0109 end -0110 if strcmpi(paramType{j},'unc') -0111 if isfield(model.annotation,'defaultLB') -0112 lb = model.annotation.defaultLB; -0113 else -0114 lb = -1000; -0115 end -0116 if isfield(model.annotation,'defaultUB') -0117 ub = model.annotation.defaultUB; -0118 else -0119 ub = 1000; -0120 end -0121 model.lb(indexes(j)) = lb; -0122 model.ub(indexes(j)) = ub; -0123 end -0124 end -0125 end -0126 if any(ismember(paramType,{'lb','ub','unc'})) -0127 invalidBound = model.lb(indexes) > model.ub(indexes); -0128 if any(invalidBound) -0129 error(['Invalid set of bounds for reaction(s): ', strjoin(model.rxns(indexes(invalidBound)),', '), '.']) -0130 end -0131 end +0097 if strcmpi(paramType{j},'var') +0098 if params(j) < 0 +0099 model.lb(indexes(j)) = params(j) * (1+var/200); +0100 model.ub(indexes(j)) = params(j) * (1-var/200); +0101 else +0102 model.lb(indexes(j)) = params(j) * (1-var/200); +0103 model.ub(indexes(j)) = params(j) * (1+var/200); +0104 end +0105 end +0106 if strcmpi(paramType{j},'unc') +0107 if isfield(model.annotation,'defaultLB') +0108 lb = model.annotation.defaultLB; +0109 else +0110 lb = -1000; +0111 end +0112 if isfield(model.annotation,'defaultUB') +0113 ub = model.annotation.defaultUB; +0114 else +0115 ub = 1000; +0116 end +0117 model.lb(indexes(j)) = lb; +0118 model.ub(indexes(j)) = ub; +0119 end +0120 end +0121 end +0122 if any(ismember(paramType,{'lb','ub','unc'})) +0123 invalidBound = model.lb(indexes) > model.ub(indexes); +0124 if any(invalidBound) +0125 error(['Invalid set of bounds for reaction(s): ', strjoin(model.rxns(indexes(invalidBound)),', '), '.']) +0126 end +0127 end
Generated by m2html © 2005
\ No newline at end of file diff --git a/doc/solver/optimizeProb.html b/doc/solver/optimizeProb.html index ebc32e00..6a0138bd 100644 --- a/doc/solver/optimizeProb.html +++ b/doc/solver/optimizeProb.html @@ -144,191 +144,194 @@

SOURCE CODE ^structUpdate(solverparams,params); -0094 -0095 % Restructering problem according to gurobi format -0096 if isfield(prob, 'csense') -0097 prob.sense = renameparams(prob.csense, {'L','G','E'}, {'<','>','='}); -0098 prob = rmfield(prob, {'csense'}); -0099 end -0100 if isfield(prob, 'osense') -0101 osense = prob.osense; -0102 prob.modelsense = renameparams(num2str(prob.osense), {'1','-1'}, {'min','max'}); -0103 prob = rmfield(prob, {'osense'}); -0104 end -0105 [prob.obj, prob.rhs, prob.vtype] = deal(prob.c, prob.b, prob.vartype); -0106 prob = rmfield(prob, {'c','b','vartype'}); -0107 -0108 resG = gurobi(prob,solverparams); -0109 -0110 try -0111 % Name output fields the same as COBRA does -0112 res.full = resG.x; -0113 res.obj = resG.objval; -0114 res.origStat = resG.status; -0115 if isfield(resG,{'pi','rc'}) -0116 res.dual = -resG.pi*osense; -0117 res.rcost = -resG.rc*osense; -0118 end -0119 if milp && strcmp(resG.status, 'TIME_LIMIT') -0120 % If res has the objval field, it succeeded, regardless of -0121 % time_limit status -0122 resG.status = 'OPTIMAL'; -0123 end -0124 switch resG.status -0125 case 'OPTIMAL' -0126 res.stat = 1; -0127 case 'UNBOUNDED' -0128 res.stat = 2; -0129 otherwise -0130 res.stat = 0; -0131 end -0132 if ~milp -0133 res.vbasis = resG.vbasis; -0134 res.cbasis = resG.cbasis; -0135 else -0136 res.mipgap = resG.mipgap; -0137 end -0138 catch -0139 res.stat = 0; -0140 res.origStat = resG.status; % useful information to have -0141 end -0142 %% Use GLPK using RAVEN-provided binary -0143 case 'glpk' -0144 solverparams.scale = 1; % Auto scaling -0145 %solverparams.tmlim = defaultparams.timeLimit; -0146 solverparams.tolbnd = defaultparams.feasTol; -0147 solverparams.toldj = defaultparams.optTol; -0148 solverparams.tolint = defaultparams.intTol; -0149 solverparams.tolobj = defaultparams.objTol; -0150 solverparams.msglev = 0; % Level of verbosity -0151 solverparams = structUpdate(solverparams,params); -0152 -0153 prob.csense = renameparams(prob.csense, {'L','G','E'}, {'U','L','S'}); -0154 -0155 if milp -0156 solverparams.tmlim = solverparams.tmlim*10; -0157 solverparams.msglev = 1; % Level of verbosity -0158 disp('Issues have been observed when using GLPK for MILP solving. Be advised to carefully observe the results, or us another solver.') -0159 end -0160 -0161 % Ensure that RAVEN glpk binary is used, return to original -0162 % directory afterwards -0163 [ravenDir,currDir]=findRAVENroot(); -0164 cd(fullfile(ravenDir,'software','GLPKmex')) -0165 [xopt, fmin, errnum, extra] = glpk(prob.c, prob.A, prob.b, prob.lb, prob.ub, prob.csense, prob.vartype, prob.osense, solverparams); -0166 cd(currDir) -0167 -0168 switch errnum % 1 = undefined; 2 = feasible; 3 = infeasible; 4 = no feasible solution; 5 = optimal; 6 = no unbounded solution -0169 case 5 -0170 res.stat = 1; % Optimal -0171 case 2 -0172 res.stat = 2; % Feasible, but not optimal -0173 otherwise -0174 res.stat = 0; -0175 end -0176 res.origStat = errnum; -0177 res.full = xopt; -0178 res.obj = fmin; -0179 res.dual = -extra.lambda*prob.osense; -0180 res.rcost = -extra.redcosts*prob.osense; -0181 %% Use scip -0182 case {'soplex','scip'} % Old 'soplex' option also allowed -0183 [xopt,fval,exitflag] = scip([], prob.c, prob.A,-prob.b, prob.b, prob.lb, prob.ub, prob.vartype); -0184 -0185 % [x,fval,exitflag,stats] = scip(H, f, A, rl, ru, lb, ub, xtype, sos, qc, nl, x0, opts) -0186 % -0187 % Input arguments*: -0188 % H - quadratic objective matrix (sparse, optional [NOT TRIL / TRIU]) -0189 % f - linear objective vector -0190 % A - linear constraint matrix (sparse) -0191 % rl - linear constraint lhs -0192 % ru - linear constraint rhs -0193 % lb - decision variable lower bounds -0194 % ub - decision variable upper bounds -0195 % xtype - string of variable integrality ('c' continuous, 'i' integer, 'b' binary) -0196 % sos - SOS structure with fields type, index and weight (see below) -0197 % qc - Quadratic Constraints structure with fields Q, l, qrl and qru (see below) -0198 % nl - Nonlinear Objective and Constraints structure (see below) -0199 % x0 - primal solution -0200 % opts - solver options (see below) -0201 % -0202 % Return arguments: -0203 % x - solution vector -0204 % fval - objective value at the solution -0205 % exitflag - exit status (see below) -0206 % stats - statistics structure -0207 % -0208 % Option Fields (all optional, see also optiset for a list): -0209 % solverOpts - specific SCIP options (list of pairs of parameter names and values) -0210 % maxiter - maximum LP solver iterations -0211 % maxnodes - maximum nodes to explore -0212 % maxtime - maximum execution time [s] -0213 % tolrfun - primal feasibility tolerance -0214 % display - solver display level [0-5] -0215 % probfile - write problem to given file -0216 % presolvedfile - write presolved problem to file -0217 % -0218 % Return Status: -0219 % 0 - Unknown -0220 % 1 - User Interrupted -0221 % 2 - Node Limit Reached -0222 % 3 - Total Node Limit Reached -0223 % 4 - Stall Node Limit Reached -0224 % 5 - Time Limit Reached -0225 % 6 - Memory Limit Reached -0226 % 7 - Gap Limit Reached -0227 % 8 - Solution Limit Reached -0228 % 9 - Solution Improvement Limit Reached -0229 % 10 - Restart Limit Reached -0230 % 11 - Problem Solved to Optimality -0231 % 12 - Problem is Infeasible -0232 % 13 - Problem is Unbounded -0233 % 14 - Problem is Either Infeasible or Unbounded -0234 -0235 res.origStat = exitflag; -0236 res.full = xopt; -0237 res.obj = fval; -0238 -0239 switch exitflag -0240 case 11 -0241 res.stat = 1; -0242 case [5, 6, 7, 8, 9, 10, 13] -0243 res.stat = 2; -0244 otherwise -0245 res.stat = 0; -0246 end -0247 otherwise -0248 error('RAVEN solver not defined or unknown. Try using setRavenSolver(''solver'').'); -0249 end -0250 if res.stat>0 -0251 res.full=res.full(1:size(prob.a,2)); +0093 if ~isempty(getCurrentTask) % If run in parallel, then one thread per gurobi +0094 solverparams.Threads=1; +0095 end +0096 solverparams = structUpdate(solverparams,params); +0097 +0098 % Restructering problem according to gurobi format +0099 if isfield(prob, 'csense') +0100 prob.sense = renameparams(prob.csense, {'L','G','E'}, {'<','>','='}); +0101 prob = rmfield(prob, {'csense'}); +0102 end +0103 if isfield(prob, 'osense') +0104 osense = prob.osense; +0105 prob.modelsense = renameparams(num2str(prob.osense), {'1','-1'}, {'min','max'}); +0106 prob = rmfield(prob, {'osense'}); +0107 end +0108 [prob.obj, prob.rhs, prob.vtype] = deal(prob.c, prob.b, prob.vartype); +0109 prob = rmfield(prob, {'c','b','vartype'}); +0110 +0111 resG = gurobi(prob,solverparams); +0112 +0113 try +0114 % Name output fields the same as COBRA does +0115 res.full = resG.x; +0116 res.obj = resG.objval; +0117 res.origStat = resG.status; +0118 if isfield(resG,{'pi','rc'}) +0119 res.dual = -resG.pi*osense; +0120 res.rcost = -resG.rc*osense; +0121 end +0122 if milp && strcmp(resG.status, 'TIME_LIMIT') +0123 % If res has the objval field, it succeeded, regardless of +0124 % time_limit status +0125 resG.status = 'OPTIMAL'; +0126 end +0127 switch resG.status +0128 case 'OPTIMAL' +0129 res.stat = 1; +0130 case 'UNBOUNDED' +0131 res.stat = 2; +0132 otherwise +0133 res.stat = 0; +0134 end +0135 if ~milp +0136 res.vbasis = resG.vbasis; +0137 res.cbasis = resG.cbasis; +0138 else +0139 res.mipgap = resG.mipgap; +0140 end +0141 catch +0142 res.stat = 0; +0143 res.origStat = resG.status; % useful information to have +0144 end +0145 %% Use GLPK using RAVEN-provided binary +0146 case 'glpk' +0147 solverparams.scale = 1; % Auto scaling +0148 %solverparams.tmlim = defaultparams.timeLimit; +0149 solverparams.tolbnd = defaultparams.feasTol; +0150 solverparams.toldj = defaultparams.optTol; +0151 solverparams.tolint = defaultparams.intTol; +0152 solverparams.tolobj = defaultparams.objTol; +0153 solverparams.msglev = 0; % Level of verbosity +0154 solverparams = structUpdate(solverparams,params); +0155 +0156 prob.csense = renameparams(prob.csense, {'L','G','E'}, {'U','L','S'}); +0157 +0158 if milp +0159 solverparams.tmlim = solverparams.tmlim*10; +0160 solverparams.msglev = 1; % Level of verbosity +0161 disp('Issues have been observed when using GLPK for MILP solving. Be advised to carefully observe the results, or us another solver.') +0162 end +0163 +0164 % Ensure that RAVEN glpk binary is used, return to original +0165 % directory afterwards +0166 [ravenDir,currDir]=findRAVENroot(); +0167 cd(fullfile(ravenDir,'software','GLPKmex')) +0168 [xopt, fmin, errnum, extra] = glpk(prob.c, prob.A, prob.b, prob.lb, prob.ub, prob.csense, prob.vartype, prob.osense, solverparams); +0169 cd(currDir) +0170 +0171 switch errnum % 1 = undefined; 2 = feasible; 3 = infeasible; 4 = no feasible solution; 5 = optimal; 6 = no unbounded solution +0172 case 5 +0173 res.stat = 1; % Optimal +0174 case 2 +0175 res.stat = 2; % Feasible, but not optimal +0176 otherwise +0177 res.stat = 0; +0178 end +0179 res.origStat = errnum; +0180 res.full = xopt; +0181 res.obj = fmin; +0182 res.dual = -extra.lambda*prob.osense; +0183 res.rcost = -extra.redcosts*prob.osense; +0184 %% Use scip +0185 case {'soplex','scip'} % Old 'soplex' option also allowed +0186 [xopt,fval,exitflag] = scip([], prob.c, prob.A,-prob.b, prob.b, prob.lb, prob.ub, prob.vartype); +0187 +0188 % [x,fval,exitflag,stats] = scip(H, f, A, rl, ru, lb, ub, xtype, sos, qc, nl, x0, opts) +0189 % +0190 % Input arguments*: +0191 % H - quadratic objective matrix (sparse, optional [NOT TRIL / TRIU]) +0192 % f - linear objective vector +0193 % A - linear constraint matrix (sparse) +0194 % rl - linear constraint lhs +0195 % ru - linear constraint rhs +0196 % lb - decision variable lower bounds +0197 % ub - decision variable upper bounds +0198 % xtype - string of variable integrality ('c' continuous, 'i' integer, 'b' binary) +0199 % sos - SOS structure with fields type, index and weight (see below) +0200 % qc - Quadratic Constraints structure with fields Q, l, qrl and qru (see below) +0201 % nl - Nonlinear Objective and Constraints structure (see below) +0202 % x0 - primal solution +0203 % opts - solver options (see below) +0204 % +0205 % Return arguments: +0206 % x - solution vector +0207 % fval - objective value at the solution +0208 % exitflag - exit status (see below) +0209 % stats - statistics structure +0210 % +0211 % Option Fields (all optional, see also optiset for a list): +0212 % solverOpts - specific SCIP options (list of pairs of parameter names and values) +0213 % maxiter - maximum LP solver iterations +0214 % maxnodes - maximum nodes to explore +0215 % maxtime - maximum execution time [s] +0216 % tolrfun - primal feasibility tolerance +0217 % display - solver display level [0-5] +0218 % probfile - write problem to given file +0219 % presolvedfile - write presolved problem to file +0220 % +0221 % Return Status: +0222 % 0 - Unknown +0223 % 1 - User Interrupted +0224 % 2 - Node Limit Reached +0225 % 3 - Total Node Limit Reached +0226 % 4 - Stall Node Limit Reached +0227 % 5 - Time Limit Reached +0228 % 6 - Memory Limit Reached +0229 % 7 - Gap Limit Reached +0230 % 8 - Solution Limit Reached +0231 % 9 - Solution Improvement Limit Reached +0232 % 10 - Restart Limit Reached +0233 % 11 - Problem Solved to Optimality +0234 % 12 - Problem is Infeasible +0235 % 13 - Problem is Unbounded +0236 % 14 - Problem is Either Infeasible or Unbounded +0237 +0238 res.origStat = exitflag; +0239 res.full = xopt; +0240 res.obj = fval; +0241 +0242 switch exitflag +0243 case 11 +0244 res.stat = 1; +0245 case [5, 6, 7, 8, 9, 10, 13] +0246 res.stat = 2; +0247 otherwise +0248 res.stat = 0; +0249 end +0250 otherwise +0251 error('RAVEN solver not defined or unknown. Try using setRavenSolver(''solver'').'); 0252 end -0253 end -0254 -0255 function s_merged=structUpdate(s_old,s_new) -0256 %Remove overlapping fields from first struct; -0257 %Obtain all unique names of remaining fields; -0258 %Merge both structs -0259 s_merged = rmfield(s_old, intersect(fieldnames(s_old), fieldnames(s_new))); -0260 names = [fieldnames(s_merged); fieldnames(s_new)]; -0261 s_merged = cell2struct([struct2cell(s_merged); struct2cell(s_new)], names, 1); -0262 end -0263 -0264 function paramlist = renameparams(paramlist,old,new) -0265 if ~iscell(paramlist) -0266 wasNoCell = true; -0267 paramlist={paramlist}; -0268 else -0269 wasNoCell = false; -0270 end -0271 for i=1:numel(old) -0272 paramlist = regexprep(paramlist,old{i},new{i}); +0253 if res.stat>0 +0254 res.full=res.full(1:size(prob.a,2)); +0255 end +0256 end +0257 +0258 function s_merged=structUpdate(s_old,s_new) +0259 %Remove overlapping fields from first struct; +0260 %Obtain all unique names of remaining fields; +0261 %Merge both structs +0262 s_merged = rmfield(s_old, intersect(fieldnames(s_old), fieldnames(s_new))); +0263 names = [fieldnames(s_merged); fieldnames(s_new)]; +0264 s_merged = cell2struct([struct2cell(s_merged); struct2cell(s_new)], names, 1); +0265 end +0266 +0267 function paramlist = renameparams(paramlist,old,new) +0268 if ~iscell(paramlist) +0269 wasNoCell = true; +0270 paramlist={paramlist}; +0271 else +0272 wasNoCell = false; 0273 end -0274 if wasNoCell -0275 paramlist=paramlist{1}; +0274 for i=1:numel(old) +0275 paramlist = regexprep(paramlist,old{i},new{i}); 0276 end -0277 end +0277 if wasNoCell +0278 paramlist=paramlist{1}; +0279 end +0280 end
Generated by m2html © 2005
\ No newline at end of file diff --git a/doc/solver/solveLP.html b/doc/solver/solveLP.html index 8c462e72..9cec93be 100644 --- a/doc/solver/solveLP.html +++ b/doc/solver/solveLP.html @@ -151,8 +151,8 @@

SOURCE CODE ^'Invalid defintion of objective function in model.c.') 0071 end 0072 %Check for valid S-matrix -0073 invalidS = ~isfinite(model.S); -0074 if any(any(invalidS)) +0073 if ~allfinite(model.S) +0074 invalidS = ~isfinite(model.S); 0075 error(['Invalid coefficients defined for reaction(s): ', strjoin(model.rxns(any(invalidS)),', '), '.']) 0076 end 0077 diff --git a/solver/optimizeProb.m b/solver/optimizeProb.m index cf33d92d..3aa3bce2 100755 --- a/solver/optimizeProb.m +++ b/solver/optimizeProb.m @@ -90,6 +90,9 @@ solverparams.FeasibilityTol = defaultparams.feasTol; solverparams.OptimalityTol = defaultparams.optTol; solverparams.Presolve = 2; + if ~isempty(getCurrentTask) % If run in parallel, then one thread per gurobi + solverparams.Threads=1; + end solverparams = structUpdate(solverparams,params); % Restructering problem according to gurobi format diff --git a/solver/solveLP.m b/solver/solveLP.m index ddd123d7..2fba3805 100755 --- a/solver/solveLP.m +++ b/solver/solveLP.m @@ -70,8 +70,8 @@ error('Invalid defintion of objective function in model.c.') end %Check for valid S-matrix -invalidS = ~isfinite(model.S); -if any(any(invalidS)) +if ~allfinite(model.S) + invalidS = ~isfinite(model.S); error(['Invalid coefficients defined for reaction(s): ', strjoin(model.rxns(any(invalidS)),', '), '.']) end