Annotation of /trunk/gadget/ecosystem.cc

Parent Directory | Revision Log

---

Revision 15 - (view) (download)

1 :	agomez	1	#include "ecosystem.h"
2 :			#include "runid.h"
3 :			#include "global.h"
4 :	ulcessvp	11	#ifndef NO_OPENMP
5 :			#include <omp.h>
6 :			#endif
7 :	agomez	1
8 :	ulcessvp	11	#ifndef NO_OPENMP
9 :			extern Ecosystem** EcoSystems;
10 :			//extern StochasticData* data;
11 :			#endif
12 :
13 :	agomez	1	Ecosystem::Ecosystem(const MainInfo& main) : printinfo(main.getPI()) {
14 :
15 :			funceval = 0;
16 :			interrupted = 0;
17 :			likelihood = 0.0;
18 :			keeper = new Keeper;
19 :
20 :			// initialise counter used when printing output files
21 :			printcount = printinfo.getPrintIteration() - 1;
22 :
23 :			// read the model specification from the main file
24 :			char* filename = main.getMainGadgetFile();
25 :			ifstream infile;
26 :			infile.open(filename, ios::in);
27 :			CommentStream commin(infile);
28 :			handle.checkIfFailure(infile, filename);
29 :			handle.Open(filename);
30 :			this->readMain(commin, main);
31 :			handle.Close();
32 :			infile.close();
33 :			infile.clear();
34 :
35 :			// if this is an optimising run then read the optimisation parameters from file
36 :			if (main.runOptimise()) {
37 :			handle.logMessage(LOGMESSAGE, ""); //write blank line to log file
38 :			if (main.getOptInfoGiven()) {
39 :			filename = main.getOptInfoFile();
40 :			infile.open(filename, ios::in);
41 :			handle.checkIfFailure(infile, filename);
42 :			handle.Open(filename);
43 :	ulcessvp	11	unsigned seed;
44 :			main.getForcePrint();
45 :			this->readOptimisation(commin, main.getSeed());
46 :	agomez	1	handle.Close();
47 :			infile.close();
48 :			infile.clear();
49 :			} else {
50 :			handle.logMessage(LOGINFO, "Warning - no optimisation file specified, using default values");
51 :			optvec.resize(new OptInfoHooke());
52 :			}
53 :			}
54 :
55 :			if (main.runOptimise())
56 :			handle.logMessage(LOGINFO, "\nFinished reading model data files, starting to run optimisation");
57 :			else
58 :			handle.logMessage(LOGINFO, "\nFinished reading model data files, starting to run simulation");
59 :			handle.logMessage(LOGMESSAGE, ""); //write blank line to log file
60 :			}
61 :
62 :			Ecosystem::~Ecosystem() {
63 :			int i;
64 :			for (i = 0; i < optvec.Size(); i++)
65 :			delete optvec[i];
66 :			for (i = 0; i < printvec.Size(); i++)
67 :			delete printvec[i];
68 :			for (i = 0; i < likevec.Size(); i++)
69 :			delete likevec[i];
70 :			for (i = 0; i < tagvec.Size(); i++)
71 :			delete tagvec[i];
72 :			for (i = 0; i < basevec.Size(); i++)
73 :			delete basevec[i];
74 :
75 :			delete Area;
76 :			delete TimeInfo;
77 :			delete keeper;
78 :			}
79 :
80 :			void Ecosystem::writeStatus(const char* filename) const {
81 :			ofstream outfile;
82 :			outfile.open(filename, ios::out);
83 :			handle.checkIfFailure(outfile, filename);
84 :			handle.Open(filename);
85 :			RUNID.Print(outfile);
86 :			outfile << "The current simulation time is " << TimeInfo->getYear()
87 :			<< ", step " << TimeInfo->getStep() << endl;
88 :
89 :			int i;
90 :			for (i = 0; i < basevec.Size(); i++)
91 :			basevec[i]->Print(outfile);
92 :			for (i = 0; i < likevec.Size(); i++)
93 :			likevec[i]->Print(outfile);
94 :			for (i = 0; i < tagvec.Size(); i++)
95 :			tagvec[i]->Print(outfile);
96 :
97 :			handle.Close();
98 :			outfile.close();
99 :			outfile.clear();
100 :			}
101 :
102 :			void Ecosystem::Reset() {
103 :			int i;
104 :			TimeInfo->Reset();
105 :			for (i = 0; i < likevec.Size(); i++)
106 :			likevec[i]->Reset(keeper);
107 :			for (i = 0; i < basevec.Size(); i++)
108 :			basevec[i]->Reset(TimeInfo);
109 :			}
110 :
111 :			void Ecosystem::Optimise() {
112 :			int i;
113 :			for (i = 0; i < optvec.Size(); i++) {
114 :	ulcessvp	11	#ifndef NO_OPENMP
115 :			int j, numThr = omp_get_max_threads ( );
116 :			DoubleVector v = this->getValues();
117 :			for (j = 0; j < numThr; j++)
118 :			EcoSystems[j]->Update(v);
119 :			#endif
120 :	ulcessvp	15	#ifdef SPECULATIVE
121 :	ulcessvp	12	optvec[i]->OptimiseLikelihoodOMP();
122 :			#else
123 :			optvec[i]->OptimiseLikelihood();
124 :			#endif
125 :	ulcessvp	11	this->writeOptValues();
126 :			}
127 :			}
128 :
129 :	agomez	1	double Ecosystem::SimulateAndUpdate(const DoubleVector& x) {
130 :			int i, j;
131 :
132 :			if (funceval == 0) {
133 :			// JMB - only need to create these vectors once
134 :			initialval.resize(keeper->numVariables(), 0.0);
135 :			currentval.resize(keeper->numVariables(), 0.0);
136 :			optflag.resize(keeper->numVariables(), 0);
137 :			keeper->getOptFlags(optflag);
138 :			}
139 :
140 :			j = 0;
141 :			keeper->getCurrentValues(currentval);
142 :			keeper->getInitialValues(initialval);
143 :			for (i = 0; i < currentval.Size(); i++) {
144 :			if (optflag[i]) {
145 :			currentval[i] = x[j] * initialval[i];
146 :			j++;
147 :			}
148 :			}
149 :
150 :			keeper->Update(currentval);
151 :			this->Simulate(0); //dont print whilst optimising
152 :
153 :			if (printinfo.getPrint()) {
154 :			printcount++;
155 :			if (printcount == printinfo.getPrintIteration()) {
156 :			keeper->writeValues(likevec, printinfo.getPrecision());
157 :			printcount = 0;
158 :			}
159 :			}
160 :
161 :			funceval++;
162 :			return likelihood;
163 :			}
164 :
165 :			void Ecosystem::writeOptValues() {
166 :			int i;
167 :			DoubleVector tmpvec(likevec.Size(), 0.0);
168 :			for (i = 0; i < likevec.Size(); i++)
169 :			tmpvec[i] = likevec[i]->getUnweightedLikelihood();
170 :
171 :	ulcessvp	11	int iters = 0;
172 :			#ifndef NO_OPENMP
173 :			int numThr = omp_get_max_threads ( );
174 :			for (i = 0; i < numThr; i++)
175 :			iters += EcoSystems[i]->getFuncEval();
176 :			#endif
177 :			handle.logMessage(LOGINFO, "\nAfter a total of", funceval+iters, "function evaluations the best point found is");
178 :	agomez	1	keeper->writeBestValues();
179 :			handle.logMessage(LOGINFO, "\nThe scores from each likelihood component are");
180 :			handle.logMessage(LOGINFO, tmpvec);
181 :			if (!isZero(keeper->getBestLikelihoodScore())) // no better point has been found
182 :			handle.logMessage(LOGINFO, "\nThe overall likelihood score is", keeper->getBestLikelihoodScore());
183 :			else
184 :			handle.logMessage(LOGINFO, "\nThe overall likelihood score is", this->getLikelihood());
185 :			}
186 :
187 :			void Ecosystem::writeInitialInformation(const char* const filename) {
188 :			keeper->openPrintFile(filename);
189 :			keeper->writeInitialInformation(likevec);
190 :			}
191 :
192 :			void Ecosystem::writeValues() {
193 :			keeper->writeValues(likevec, printinfo.getPrecision());
194 :			}
195 :
196 :			void Ecosystem::writeParams(const char* const filename, int prec) const {
197 :			if ((funceval > 0) && (interrupted == 0)) {
198 :			//JMB - print the final values to any output files specified
199 :			//in case they have been missed by the -print value
200 :			if (printinfo.getPrint())
201 :			keeper->writeValues(likevec, printinfo.getPrecision());
202 :			}
203 :			keeper->writeParams(optvec, filename, prec, interrupted);
204 :			}

---

root@forge.cesga.es	ViewVC Help
Powered by ViewVC 1.0.0