defineLinearConstraints

PURPOSE ^

SYNOPSIS ^

function [A, b_L, b_U, model] = defineLinearConstraints(model, method)

DESCRIPTION ^

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SOURCE CODE ^

0001 function [A, b_L, b_U, model] = defineLinearConstraints(model, method)
0002 
0003 fOffset = 1e-4;
0004 
0005 if nargin < 2
0006     method = 1;
0007 end
0008 
0009 b_L = model.lb;
0010 b_U = model.ub;
0011 A = model.N;
0012 
0013 %marked = zeros(length(model.lb), 1);
0014 % for i = 1:size(A, 1)
0015 %     if all(abs(A(i,:)) < 1e-5)
0016 %         marked(i) = 2;
0017 %         if b_L(i) > 1e-5 || b_U(i) < -1e-5
0018 %             display('upperbound/lowerbound error');
0019 %             pause;
0020 %         end
0021 %     end
0022 % end
0023 % marked
0024 % pause;
0025 
0026 %marked = zeros(size(A,1), 1);
0027 % for  i = 1:length(model.lb)
0028 %     if mod(i,40) == 0, i, end
0029 %
0030 %     tmarked = marked;
0031 %     tmarked(i) = 1;
0032 %
0033 %     LPproblem.A = model.S;
0034 %
0035 %     LPproblem.b = zeros(length(model.mets),1);
0036 %     LPproblem.csense = ['E' * ones(length(model.mets), 1)];
0037 %     LPproblem.ub = model.ub;
0038 %     LPproblem.ub(tmarked ~= 0) = 1e6;
0039 %     LPproblem.lb = model.lb;
0040 %     LPproblem.lb(tmarked ~= 0) = -1e6;
0041 %
0042 %     LPproblem.osense = 1;
0043 %
0044 %     %minimize
0045 %     LPproblem.c = zeros(length(model.lb),1);
0046 %     LPproblem.c(i) = 1;
0047 %     soln = solveCobraLP(LPproblem);
0048 %     vmin = soln.obj;
0049 %
0050 %     %maximize
0051 %     LPproblem.c = -LPproblem.c;
0052 %     soln = solveCobraLP(LPproblem);
0053 %     vmax = -soln.obj;
0054 %
0055 %     %[i b_L(i) vmin vmax b_U(i)]
0056 %     if (vmin > b_L(i) + 1e-2 && vmax < b_U(i) - 1e-2)
0057 %         %[i b_L(i) vmin vmax b_U(i)]
0058 %         marked(i) = 1;
0059 %     end
0060 %
0061 % end
0062 % display ('done1');
0063 
0064 marked = zeros(size(model.lb));
0065 n = length(model.lb);
0066 dirs = zeros(size(model.lb));
0067 for i = 1:n
0068     if model.lb(i) == 0 && model.ub(i) > fOffset
0069         dirs(i)  = 1;
0070     elseif model.lb(i) < -fOffset && model.ub(i) == 0
0071         dirs(i) = -1;
0072     end
0073 end
0074 
0075 
0076 
0077 for i = 1:n
0078     if mod(i,10) == 0
0079         i;
0080     end
0081     
0082     LPproblem.A = model.S;
0083 
0084     LPproblem.b = zeros(length(model.mets),1);
0085     LPproblem.csense = ['E'*ones(length(model.mets), 1)];
0086     LPproblem.csense = char(LPproblem.csense);
0087     LPproblem.ub = model.ub;
0088     LPproblem.lb = model.lb;
0089     LPproblem.lb(marked > 0) = fOffset;
0090     LPproblem.ub(marked < 0) = -fOffset;
0091     
0092     LPproblem.osense = 1;
0093 
0094     %maximize
0095     LPproblem.c = -dirs;
0096     LPproblem.c(marked == 0 ) = 0;
0097     soln = solveCobraLP(LPproblem);
0098     vm = soln.full;
0099     
0100     if soln.stat ~= 1
0101         soln
0102         i
0103         pause;
0104     end
0105     
0106     [model.lb model.ub dirs vm];
0107 
0108     indexes = (vm.*dirs > fOffset) & marked == 0 & dirs ~= 0;
0109     markednew = marked;
0110     markednew(indexes) = dirs(indexes);
0111     
0112     [model.lb model.ub dirs vm markednew];
0113 
0114     if markednew == marked
0115         break;
0116     else
0117         marked = markednew;
0118     end
0119 %     if( vmin < b_L(i) - 1e-7 || vmax > b_U(i) + 1e-7)
0120 %         [i b_L(i) vmin vmax b_U(i)]
0121 %         pause;
0122 %     end
0123 %     if (abs(vmin) < 1e-6 && abs(vmax) < 1e-6)
0124 %         marked(i) = 2;
0125 %     end
0126 end
0127 %model.redundant_constraints = marked;
0128 if method == 1
0129     if isfield(model, 'N')
0130         A = model.N;
0131     else
0132         A = null(model.S);
0133     end
0134 elseif method ==2 
0135     A = model.S;
0136 end
0137 b_L(marked > 0) = fOffset;
0138 b_U(marked < 0) = -fOffset;
0139 [b_L b_U];
0140 x = sum(abs(A),2) < 1e-9;
0141 b_L = b_L(~x);
0142 b_U = b_U(~x);
0143 A = A(~x, :);
0144 return

Generated on Thu 21-Jun-2012 15:39:23 by m2html © 2003