Main Page   Namespace List   Class Hierarchy   Alphabetical List   Compound List   File List   Compound Members   File Members  

EGFRTraffickingNIH3T3MinimizableDirector.cpp

Go to the documentation of this file.
00001 // EGFRTraffickingNIH3T3MinimizableDirector.cpp: implementation of the CEGFRTraffickingNIH3T3MinimizableDirector class.
00002 //
00004 
00005 #include "EGFRTraffickingNIH3T3MinimizableDirector.h"
00006 
00008 // Construction/Destruction
00010 
00011 CEGFRTraffickingNIH3T3MinimizableDirector::CEGFRTraffickingNIH3T3MinimizableDirector()
00012 {
00013         // create networks,experiments,observers,movers,minimizables
00014         double maxSimTime = 361.0;
00015         double stepSize = 0.0001;
00016         double frequency = 1.0;
00017         int nCells = 1;
00018         double moveTol = 1.0e-06;
00019 
00020         m_dGammaSquared = 0.0;
00021         
00022         // flags for all minimizables
00023         bool rateFlag = true;
00024         bool chemFlag = false;
00025         bool logsFlag = false;
00026         double timeWeight = 1.0;
00027         double rateWeight = 0.0;
00028         double chemWeight = 0.0;
00029 
00030         int cnt = 0;
00031         /*
00032         // NIH-3T3 cells
00033         // expt. 0 - NIH-3T3 cells transfected with 1.0(5) X 10^5 receptors, treated 
00034         //              with EGF at 100 nM concentration
00035         _networkList.push_back(new CEGFRTraffickingNetwork());
00036         _experimentList.push_back(new Experiment(_networkList[cnt],"EGF at 100 nM, Helin 1991"));
00037         // put data files here
00038         _experimentList[cnt]->AddDataFileToList("../DataFiles/TimeSeries/SurfaceEGFREGFHelin100nM.dat");
00039         _experimentList[cnt]->AddForcingFileToList("../DataFiles/TimeSeries/TwoHourEGF.dat");
00040         // mover and plotter
00041         //_moverList.push_back(new CQualityControlCashKarpMover(frequency,stepSize,moveTol));
00042         _moverList.push_back(new CStiffBulirschStoerMover(frequency,stepSize,moveTol));
00043         // minimizable
00044         _minimizableList.push_back(new SingleNetworkMinimizable());
00045         // parameters
00046         _minimizableList[cnt]->Initialize(_experimentList[cnt],_moverList[cnt],nCells,rateFlag,chemFlag,logsFlag,timeWeight,rateWeight,chemWeight);
00047         _minimizableList[cnt]->GetConversionFactor(67)->SetFactorFixed(true);
00048         _minimizableList[cnt]->GetConversionFactor(67)->SetFactorValue(100.0/100000.0);
00049         _minimizableList[cnt]->SetRateConstantsWeight(0.0);
00050         // attach observer
00051         
00052         _plotterList.push_back(new GnuPlotterTimeSeriesObserver());
00053         _plotterList[cnt]->SetTitle("EGFR Trafficking Network:EGF at 100 nM (Helin 1991):ODE");
00054         _minimizableList[cnt]->Attach(_plotterList[cnt]);
00055         
00056         cnt++;
00057         */
00058         
00059         // expt. 1 - NIH-3T3's transfected w/ wild-type EGFR; measurements made
00060         //              with 10-30% of receptors initially occupied; cells express 
00061         //              1.0(5) X 10^5 receptors
00062         _networkList.push_back(new CEGFRTraffickingNetwork());
00063         _experimentList.push_back(new Experiment(_networkList[cnt],"20 Percent Initial Occupancy, Helin 1991"));
00064         // put data files here
00065         _experimentList[cnt]->AddDataFileToList("../DataFiles/TimeSeries/EGFHelinOccupied.dat");
00066         // mover and plotter
00067         //_moverList.push_back(new CQualityControlCashKarpMover(frequency,stepSize,moveTol));
00068         _moverList.push_back(new CStiffBulirschStoerMover(frequency,stepSize,moveTol));
00069         // minimizable
00070         _minimizableList.push_back(new SingleNetworkMinimizable());
00071         // parameters
00072         _minimizableList[cnt]->Initialize(_experimentList[cnt],_moverList[cnt],nCells,rateFlag,chemFlag,logsFlag,timeWeight,rateWeight,chemWeight);
00073         _minimizableList[cnt]->GetConversionFactor(68)->SetFactorFixed(true);
00074         _minimizableList[cnt]->GetConversionFactor(68)->SetFactorValue(100.0/20000.0);
00075         _minimizableList[cnt]->SetRateConstantsWeight(0.0);
00076         // attach observer
00077         
00078         _plotterList.push_back(new GnuPlotterTimeSeriesObserver());
00079         _plotterList[cnt]->SetTitle("EGFR Trafficking Network:20 Percent Initial EGFR Occupancy (Helin 1991):ODE");
00080         _minimizableList[cnt]->Attach(_plotterList[cnt]);
00081         
00082         cnt++;
00083 
00084         // expt. 2 - NIH-3T3's transfected with wt EGFR; 1 X 10^5 receptors, EGF of
00085         //              100 ng/ml
00086         _networkList.push_back(new CEGFRTraffickingNetwork());
00087         _experimentList.push_back(new Experiment(_networkList[cnt],"EGF at 100 ng/ml wt1, Carter 1998"));
00088         // put data files here
00089         _experimentList[cnt]->AddDataFileToList("../DataFiles/TimeSeries/SurfaceEGFREGFwt1Carter.dat");
00090         _experimentList[cnt]->AddForcingFileToList("../DataFiles/TimeSeries/TwoHourEGF.dat");
00091         // mover and plotter
00092         //_moverList.push_back(new CQualityControlCashKarpMover(frequency,stepSize,moveTol));
00093         _moverList.push_back(new CStiffBulirschStoerMover(frequency,stepSize,moveTol));
00094         // minimizable
00095         _minimizableList.push_back(new SingleNetworkMinimizable());
00096         // parameters
00097         _minimizableList[cnt]->Initialize(_experimentList[cnt],_moverList[cnt],nCells,rateFlag,chemFlag,logsFlag,timeWeight,rateWeight,chemWeight);
00098         _minimizableList[cnt]->GetConversionFactor(67)->SetFactorFixed(true);
00099         _minimizableList[cnt]->GetConversionFactor(67)->SetFactorValue(100.0/100000.0);
00100         _minimizableList[cnt]->SetRateConstantsWeight(0.0);
00101         // attach observer
00102         
00103         _plotterList.push_back(new GnuPlotterTimeSeriesObserver());
00104         _plotterList[cnt]->SetTitle("EGFR Trafficking Network:EGF at 100 ng/ml, wt1 (Carter 1998):ODE");
00105         _minimizableList[cnt]->Attach(_plotterList[cnt]);
00106         
00107         cnt++;
00108 
00109         // expt. 3 - NIH-3T3's transfected with wt EGFR; 2.75 X 10^5 receptors, EGF of
00110         //              100 ng/ml
00111         _networkList.push_back(new CEGFRTraffickingNetwork());
00112         _experimentList.push_back(new Experiment(_networkList[cnt],"EGF at 100 ng/ml wt2, Carter 1998"));
00113         // put data files here
00114         _experimentList[cnt]->AddDataFileToList("../DataFiles/TimeSeries/SurfaceEGFREGFwt2Carter.dat");
00115         _experimentList[cnt]->AddForcingFileToList("../DataFiles/TimeSeries/TwoHourEGF.dat");
00116         // mover and plotter
00117         //_moverList.push_back(new CQualityControlCashKarpMover(frequency,stepSize,moveTol));
00118         _moverList.push_back(new CStiffBulirschStoerMover(frequency,stepSize,moveTol));
00119         // minimizable
00120         _minimizableList.push_back(new SingleNetworkMinimizable());
00121         // parameters
00122         _minimizableList[cnt]->Initialize(_experimentList[cnt],_moverList[cnt],nCells,rateFlag,chemFlag,logsFlag,timeWeight,rateWeight,chemWeight);
00123         _minimizableList[cnt]->GetConversionFactor(67)->SetFactorFixed(true);
00124         _minimizableList[cnt]->GetConversionFactor(67)->SetFactorValue(100.0/275000.0);
00125         _minimizableList[cnt]->SetRateConstantsWeight(0.0);
00126         // attach observer
00127         
00128         _plotterList.push_back(new GnuPlotterTimeSeriesObserver());
00129         _plotterList[cnt]->SetTitle("EGFR Trafficking Network:EGF at 100 ng/ml, wt2 (Carter 1998):ODE");
00130         _minimizableList[cnt]->Attach(_plotterList[cnt]);
00131         
00132         cnt++;
00133         
00134         // expt. 4 - NIH-3T3's with native EGFR (Jin's data). EGF of 100 ng/ml
00135         _networkList.push_back(new CEGFRTraffickingNetwork());
00136         _experimentList.push_back(new Experiment(_networkList[cnt],"EGF at 100 ng/ml (Wu 2002)"));
00137         // put data files here
00138         _experimentList[cnt]->AddDataFileToList("../DataFiles/TimeSeries/TotalEGFRWuVector.dat");
00139         _experimentList[cnt]->AddDataFileToList("../DataFiles/TimeSeries/ErkEGFWuVector.dat");
00140         _experimentList[cnt]->AddForcingFileToList("../DataFiles/TimeSeries/TwoHourEGF.dat");
00141         // mover and plotter
00142         //_moverList.push_back(new CQualityControlCashKarpMover(frequency,stepSize,moveTol));
00143         _moverList.push_back(new CStiffBulirschStoerMover(frequency,stepSize,moveTol));
00144         // minimizable
00145         _minimizableList.push_back(new SingleNetworkMinimizable());
00146         // parameters
00147         _minimizableList[cnt]->Initialize(_experimentList[cnt],_moverList[cnt],nCells,rateFlag,chemFlag,logsFlag,timeWeight,rateWeight,chemWeight);
00148         _minimizableList[cnt]->GetConversionFactor(77)->SetFactorFixed(true);
00149         _minimizableList[cnt]->GetConversionFactor(77)->SetFactorValue(1.0/8000.0);
00150         _minimizableList[cnt]->SetRateConstantsWeight(0.0);
00151         // attach observer
00152         
00153         _plotterList.push_back(new GnuPlotterTimeSeriesObserver());
00154         _plotterList[cnt]->SetTitle("EGFR Trafficking Network:EGF at 100 ng/ml, Endogenous EGFR (Wu 2002):ODE");
00155         _minimizableList[cnt]->Attach(_plotterList[cnt]);
00156         
00157         cnt++;
00158 
00159         // expt. 4 - NIH-3T3's with native EGFR and stable F28L Cdc42 (Jin's data). 
00160         //              EGF of 100 ng/ml
00161         _networkList.push_back(new CEGFRTraffickingNetwork());
00162         _experimentList.push_back(new Experiment(_networkList[cnt],"EGF at 100 ng/ml, stable F28L (Wu 2002)"));
00163         // put data files here
00164         _experimentList[cnt]->AddDataFileToList("../DataFiles/TimeSeries/TotalEGFRWuF28L.dat");
00165         _experimentList[cnt]->AddDataFileToList("../DataFiles/TimeSeries/ErkEGFWuF28L.dat");
00166         _experimentList[cnt]->AddForcingFileToList("../DataFiles/TimeSeries/TwoHourEGF.dat");
00167         // mover and plotter
00168         //_moverList.push_back(new CQualityControlCashKarpMover(frequency,stepSize,moveTol));
00169         _moverList.push_back(new CStiffBulirschStoerMover(frequency,stepSize,moveTol));
00170         // minimizable
00171         _minimizableList.push_back(new SingleNetworkMinimizable());
00172         // parameters
00173         _minimizableList[cnt]->Initialize(_experimentList[cnt],_moverList[cnt],nCells,rateFlag,chemFlag,logsFlag,timeWeight,rateWeight,chemWeight);
00174         _minimizableList[cnt]->GetConversionFactor(77)->SetFactorFixed(true);
00175         _minimizableList[cnt]->GetConversionFactor(77)->SetFactorValue(1.0/8000.0);
00176         _minimizableList[cnt]->SetRateConstantsWeight(0.0);
00177         // attach observer
00178         
00179         _plotterList.push_back(new GnuPlotterTimeSeriesObserver());
00180         _plotterList[cnt]->SetTitle("EGFR Trafficking Network:EGF at 100 ng/ml, stable F28L (Wu 2002):ODE");
00181         _minimizableList[cnt]->Attach(_plotterList[cnt]);
00182         
00183         cnt++;
00184 
00185         // expt. 4 - NIH-3T3's with native EGFR and stable F28L/dL8 Cdc42 (Jin's data). 
00186         //              EGF of 100 ng/ml
00187         _networkList.push_back(new CEGFRTraffickingNetwork());
00188         _experimentList.push_back(new Experiment(_networkList[cnt],"EGF at 100 ng/ml, stable F28L/dL8 (Wu 2002)"));
00189         // put data files here
00190         _experimentList[cnt]->AddDataFileToList("../DataFiles/TimeSeries/TotalEGFRWuF28LdL8.dat");
00191         _experimentList[cnt]->AddDataFileToList("../DataFiles/TimeSeries/ErkEGFWuF28LdL8.dat");
00192         _experimentList[cnt]->AddForcingFileToList("../DataFiles/TimeSeries/TwoHourEGF.dat");
00193         // mover and plotter
00194         //_moverList.push_back(new CQualityControlCashKarpMover(frequency,stepSize,moveTol));
00195         _moverList.push_back(new CStiffBulirschStoerMover(frequency,stepSize,moveTol));
00196         // minimizable
00197         _minimizableList.push_back(new SingleNetworkMinimizable());
00198         // parameters
00199         _minimizableList[cnt]->Initialize(_experimentList[cnt],_moverList[cnt],nCells,rateFlag,chemFlag,logsFlag,timeWeight,rateWeight,chemWeight);
00200         _minimizableList[cnt]->GetConversionFactor(77)->SetFactorFixed(true);
00201         _minimizableList[cnt]->GetConversionFactor(77)->SetFactorValue(1.0/8000.0);
00202         _minimizableList[cnt]->SetRateConstantsWeight(0.0);
00203         // attach observer
00204         
00205         _plotterList.push_back(new GnuPlotterTimeSeriesObserver());
00206         _plotterList[cnt]->SetTitle("EGFR Trafficking Network:EGF at 100 ng/ml, stable F28L/dL8 (Wu 2002):ODE");
00207         _minimizableList[cnt]->Attach(_plotterList[cnt]);
00208         
00209         cnt++;
00210 
00211         // parameters which need a prior
00212         //m_iPriorList.push_back(1);
00213         //m_iPriorList.push_back(2);
00214         //m_iPriorList.push_back(11);
00215         //m_iPriorList.push_back(24);
00216 
00217         // allocate the master residuals list
00218         int masterResidualsSize = 0;
00219         for(int i = 0; i < _minimizableList.size(); i++)
00220         {
00221                 int subRes = _minimizableList[i]->GetNResiduals();
00222                 masterResidualsSize += subRes;
00223                 cout << subRes << " residuals in minimizable " << i << endl;
00224         }
00225         // add residuals for prior parameter costs, if applicable
00226         masterResidualsSize += m_iPriorList.size();
00227         cout << m_iPriorList.size() << " parameters are constrained by a prior " << endl;
00228         if(m_dGammaSquared > 0.0)
00229         {
00230                 masterResidualsSize += 1;
00231         }
00232         Allocate(masterResidualsSize);
00233 
00234         cout << "Master Residuals size : " << masterResidualsSize << endl;
00235 
00236         // define the experiments
00237         DefineExperiments();
00238 
00239         //DumpResidualInfo();
00240 
00241 }
00242 
00243 CEGFRTraffickingNIH3T3MinimizableDirector::~CEGFRTraffickingNIH3T3MinimizableDirector()
00244 {
00245         for(int i = 0; i < _plotterList.size(); i++)
00246         {
00247                 _minimizableList[i]->Detach(_plotterList[i]);
00248                 delete _plotterList[i];
00249         }
00250         m_iPriorList.erase(m_iPriorList.begin(),m_iPriorList.end());
00251 }
00252 
00253 void CEGFRTraffickingNIH3T3MinimizableDirector::DefineExperiments()
00254 {       
00255         /*
00256         // 3T3 cells
00257         // 100,000 receptors, EGF at 100 nM
00258         _minimizableList[0]->GetExperiment()->GetReactionNetwork()->RemoveReaction(115);
00259         _minimizableList[0]->GetExperiment()->GetReactionNetwork()->GetChemical(0)->SetInitialAmount(7.2e+007);
00260         _minimizableList[0]->GetExperiment()->GetReactionNetwork()->GetChemical(3)->SetInitialAmount(10000);
00261         _minimizableList[0]->GetExperiment()->GetReactionNetwork()->GetChemical(4)->SetInitialAmount(90000);    
00262         _minimizableList[0]->GetExperiment()->GetReactionNetwork()->ChemicalReset();
00263         */
00264         
00265         // No EGF, 20% of receptors on surface initially occupied
00266         _minimizableList[0]->GetExperiment()->GetReactionNetwork()->RemoveReaction(115);
00267         _minimizableList[0]->GetExperiment()->GetReactionNetwork()->GetChemical(3)->SetInitialAmount(10000.0);
00268         _minimizableList[0]->GetExperiment()->GetReactionNetwork()->GetChemical(4)->SetInitialAmount(70000.0);
00269         _minimizableList[0]->GetExperiment()->GetReactionNetwork()->GetChemical(5)->SetInitialAmount(10000.0);
00270         _minimizableList[0]->GetExperiment()->GetReactionNetwork()->GetChemical(6)->SetInitialAmount(10000.0);
00271         _minimizableList[0]->GetExperiment()->GetReactionNetwork()->ChemicalReset();
00272 
00273         // 100,000 receptors assumed, EGF at 100 ng/ml
00274         _minimizableList[1]->GetExperiment()->GetReactionNetwork()->RemoveReaction(115);
00275         _minimizableList[1]->GetExperiment()->GetReactionNetwork()->GetChemical(0)->SetInitialAmount(1.2e+007);
00276         _minimizableList[1]->GetExperiment()->GetReactionNetwork()->GetChemical(3)->SetInitialAmount(10000);
00277         _minimizableList[1]->GetExperiment()->GetReactionNetwork()->GetChemical(4)->SetInitialAmount(90000);
00278         _minimizableList[1]->GetExperiment()->GetReactionNetwork()->ChemicalReset();
00279 
00280         // 275,000 receptors assumed, EGF at 100 ng/ml
00281         _minimizableList[2]->GetExperiment()->GetReactionNetwork()->RemoveReaction(115);
00282         _minimizableList[2]->GetExperiment()->GetReactionNetwork()->GetChemical(0)->SetInitialAmount(1.2e+007);
00283         _minimizableList[2]->GetExperiment()->GetReactionNetwork()->GetChemical(3)->SetInitialAmount(10000);
00284         _minimizableList[2]->GetExperiment()->GetReactionNetwork()->GetChemical(4)->SetInitialAmount(265000);
00285         _minimizableList[2]->GetExperiment()->GetReactionNetwork()->ChemicalReset();
00286         
00287         // 8,000 receptors assumed, EGF at 100 ng/ml
00288         _minimizableList[3]->GetExperiment()->GetReactionNetwork()->RemoveReaction(115);
00289         _minimizableList[3]->GetExperiment()->GetReactionNetwork()->GetChemical(0)->SetInitialAmount(1.2e+007);
00290         _minimizableList[3]->GetExperiment()->GetReactionNetwork()->GetChemical(3)->SetInitialAmount(8000);
00291         _minimizableList[3]->GetExperiment()->GetReactionNetwork()->GetChemical(4)->SetInitialAmount(0.0);
00292         _minimizableList[3]->GetExperiment()->GetReactionNetwork()->ChemicalReset();
00293 
00294         // 8,000 receptors assumed, EGF at 100 ng/ml, F28L 42 5-fold overexpressed
00295         _minimizableList[4]->GetExperiment()->GetReactionNetwork()->GetChemical(0)->SetInitialAmount(0.0);
00296         _minimizableList[4]->GetExperiment()->GetReactionNetwork()->GetChemical(3)->SetInitialAmount(8000);
00297         _minimizableList[4]->GetExperiment()->GetReactionNetwork()->GetChemical(4)->SetInitialAmount(0.0);
00298         // start with overexpressed, all-active 42
00299         _minimizableList[4]->GetExperiment()->GetReactionNetwork()->GetChemical(46)->SetInitialAmount(0.0);
00300         _minimizableList[4]->GetExperiment()->GetReactionNetwork()->GetChemical(47)->SetInitialAmount(0.5*600000);
00301         _minimizableList[4]->GetExperiment()->GetReactionNetwork()->GetChemical(48)->SetInitialAmount(0.0);
00302         _minimizableList[4]->GetExperiment()->GetReactionNetwork()->GetChemical(49)->SetInitialAmount(0.5*600000);
00303         _minimizableList[4]->GetExperiment()->GetReactionNetwork()->ChemicalReset();
00304         
00305         // 8,000 receptors assumed, EGF at 100 ng/ml, F28L/dL8 42 5-fold overexpressed
00306         _minimizableList[5]->GetExperiment()->GetReactionNetwork()->RemoveReaction(53);
00307         _minimizableList[5]->GetExperiment()->GetReactionNetwork()->RemoveReaction(51);
00308         _minimizableList[5]->GetExperiment()->GetReactionNetwork()->GetChemical(0)->SetInitialAmount(0.0);
00309         _minimizableList[5]->GetExperiment()->GetReactionNetwork()->GetChemical(3)->SetInitialAmount(8000);
00310         _minimizableList[5]->GetExperiment()->GetReactionNetwork()->GetChemical(4)->SetInitialAmount(0.0);
00311         // start with all 42 active
00312         _minimizableList[5]->GetExperiment()->GetReactionNetwork()->GetChemical(46)->SetInitialAmount(0.0);
00313         _minimizableList[5]->GetExperiment()->GetReactionNetwork()->GetChemical(47)->SetInitialAmount(0.5*600000);
00314         _minimizableList[5]->GetExperiment()->GetReactionNetwork()->GetChemical(48)->SetInitialAmount(0.0);
00315         _minimizableList[5]->GetExperiment()->GetReactionNetwork()->GetChemical(49)->SetInitialAmount(0.5*600000);
00316         _minimizableList[5]->GetExperiment()->GetReactionNetwork()->ChemicalReset();
00317 }
00318 
00319 double CEGFRTraffickingNIH3T3MinimizableDirector::GetParameter(int parIndex)
00320 {
00321         // for this director, everyone has the same parameters so any set 
00322         // can be returned
00323         //return _minimizableList[0]->GetParameter(parIndex);
00324         
00325         if(parIndex < _minimizableList[0]->GetNParameters())
00326         {
00327                 return _minimizableList[0]->GetParameter(parIndex);
00328         }
00329 }
00330 
00331 int CEGFRTraffickingNIH3T3MinimizableDirector::GetNParameters()
00332 {
00333         // ditto as for GetParameters()
00334         //return _minimizableList[0]->GetNParameters();
00335 
00336         // add in chemical concentrations being optimized
00337         return _minimizableList[0]->GetNParameters() + 9;
00338 }
00339 
00340 double CEGFRTraffickingNIH3T3MinimizableDirector::ComputeResiduals(double *parameters)
00341 {
00342         int lastrc = _minimizableList[0]->GetNParameters();
00343         // set all the ICs first
00344         for(int i = 0; i < 4; i++)
00345         {
00346                 _minimizableList[i]->GetExperiment()->GetReactionNetwork()->GetChemical(45)->SetInitialAmount(parameters[lastrc]);
00347                 _minimizableList[i]->GetExperiment()->GetReactionNetwork()->GetChemical(46)->SetInitialAmount(parameters[lastrc+1]);
00348                 _minimizableList[i]->GetExperiment()->GetReactionNetwork()->GetChemical(48)->SetInitialAmount(parameters[lastrc+2]);
00349                 _minimizableList[i]->GetExperiment()->GetReactionNetwork()->GetChemical(51)->SetInitialAmount(parameters[lastrc+3]);
00350                 _minimizableList[i]->GetExperiment()->GetReactionNetwork()->GetChemical(57)->SetInitialAmount(parameters[lastrc+4]);
00351                 _minimizableList[i]->GetExperiment()->GetReactionNetwork()->GetChemical(58)->SetInitialAmount(parameters[lastrc+5]);
00352                 _minimizableList[i]->GetExperiment()->GetReactionNetwork()->GetChemical(63)->SetInitialAmount(parameters[lastrc+6]);
00353                 _minimizableList[i]->GetExperiment()->GetReactionNetwork()->GetChemical(71)->SetInitialAmount(parameters[lastrc+7]);    
00354                 _minimizableList[i]->GetExperiment()->GetReactionNetwork()->ChemicalReset();
00355 
00356         }
00357         for(i = 4; i < 6; i++)
00358         {
00359                 double mult = 1.0 + parameters[lastrc + 8];
00360                 _minimizableList[i]->GetExperiment()->GetReactionNetwork()->GetChemical(45)->SetInitialAmount(parameters[lastrc]);
00361                 _minimizableList[i]->GetExperiment()->GetReactionNetwork()->GetChemical(46)->SetInitialAmount(0.0);
00362                 _minimizableList[i]->GetExperiment()->GetReactionNetwork()->GetChemical(48)->SetInitialAmount(0.0);
00363                 _minimizableList[i]->GetExperiment()->GetReactionNetwork()->GetChemical(51)->SetInitialAmount(parameters[lastrc+3]);
00364                 _minimizableList[i]->GetExperiment()->GetReactionNetwork()->GetChemical(57)->SetInitialAmount(parameters[lastrc+4]);
00365                 _minimizableList[i]->GetExperiment()->GetReactionNetwork()->GetChemical(58)->SetInitialAmount(parameters[lastrc+5]);
00366                 _minimizableList[i]->GetExperiment()->GetReactionNetwork()->GetChemical(63)->SetInitialAmount(parameters[lastrc+6]);    
00367                 _minimizableList[i]->GetExperiment()->GetReactionNetwork()->GetChemical(71)->SetInitialAmount(parameters[lastrc+6]);    
00368                 _minimizableList[i]->GetExperiment()->GetReactionNetwork()->GetChemical(47)->SetInitialAmount(mult*parameters[lastrc+1]);
00369                 _minimizableList[i]->GetExperiment()->GetReactionNetwork()->GetChemical(49)->SetInitialAmount(mult*parameters[lastrc+2]);
00370                 _minimizableList[i]->GetExperiment()->GetReactionNetwork()->ChemicalReset();
00371         }
00372 
00373         // now fire away
00374         // just sum up the costs from each director and then add the 
00375         // integration penalty, = (gamma^2/2)*(sum_ti)^2
00376         double summedCost = 0.0;
00377         double timeCost = 0.0;
00378         double pcost = 0.0;
00379         
00380         // compute the residuals/cost from each minimizable
00381         for(int i = 0; i < _minimizableList.size(); i++)
00382         {
00383                 summedCost += _minimizableList[i]->ObjectiveFunction(parameters);
00384         }
00385         
00386         // load up the master list with each set of residuals
00387         int nextResidual = 0;
00388         for(i = 0; i < _minimizableList.size(); i++)
00389         {
00390                 int nRes = _minimizableList[i]->GetNResiduals();
00391                 for(int j = 0; j < nRes; j++)
00392                 {
00393                         residuals[nextResidual+j] = (_minimizableList[i]->GetResiduals())[j];
00394                 }
00395                 nextResidual += nRes;
00396         }
00397         
00398         // now compute the cost for deviation from guessed values, if any
00399         for(i = 0; i < m_iPriorList.size(); i++)
00400         {
00401                 int which = m_iPriorList[i];
00402                 double presid = _networkList[0]->GetRateConstant(which)->GetValue() - _networkList[0]->GetRateConstant(which)->GetInitialValue();
00403                 presid = presid/(sqrt(2.0)*_networkList[0]->GetRateConstant(which)->GetErrorInInitialValue());
00404                 residuals[nextResidual + i] = presid;
00405                 pcost += presid*presid;
00406         }
00407         nextResidual += m_iPriorList.size();
00408         
00409         // if the prefactor for the integration cost is > 0, compute its residual (otherwise
00410         // that residual shouldn't exist)
00411         if(m_dGammaSquared > 0.0)
00412         {
00413                 double totalSteps = 0.0;
00414                 for(int l = 0; l < _moverList.size(); l++)
00415                 {
00416                         totalSteps += _moverList[l]->GetTotalStepsTaken();
00417                 }
00418                 residuals[nResiduals - 1] = sqrt(m_dGammaSquared/2.0)*totalSteps;
00419                 timeCost = residuals[nResiduals - 1]*residuals[nResiduals - 1];
00420         }
00421         return summedCost + pcost + timeCost;
00422 }
Generated on Mon Nov 3 09:37:47 2003 by doxygen1.2.14 written by Dimitri van Heesch, © 1997-2002