1919TransientInput::TransientInput (){/* {{{*/
2020
2121 this ->enum_type =UNDEF;
22- this ->inputs =NULL ;
23- this ->numtimesteps =0 ;
2422 this ->parameters =NULL ;
25- this ->timesteps =NULL ;
2623
2724 this ->current_input =NULL ;
2825 this ->current_step =-1 ;
@@ -38,18 +35,8 @@ TransientInput::TransientInput(int in_enum_type,int nbe,int nbv,IssmDouble* time
3835
3936 /* Allocate values and timesteps, and copy: */
4037 _assert_ (N>=0 && N<1e6 );
41- this ->numtimesteps =N;
42- if (N>0 ){
43- this ->timesteps =xNew<IssmDouble>(N);
44- xMemCpy (this ->timesteps ,timesin,N);
45-
46- this ->inputs = xNew<Input*>(N);
47- for (int i=0 ;i<N;i++) this ->inputs [i] = NULL ;
48- }
49- else {
50- this ->timesteps =0 ;
51- this ->inputs =0 ;
52- }
38+ this ->numtimesteps = N;
39+ for (int i=0 ;i<N;i++) this ->timesteps .push_back (timesin[i]);
5340 this ->parameters = NULL ;
5441 this ->current_input =NULL ;
5542 this ->current_step =-1 ;
@@ -60,8 +47,8 @@ TransientInput::~TransientInput(){/*{{{*/
6047 for (int i=0 ;i<this ->numtimesteps ;i++){
6148 delete this ->inputs [i];
6249 }
63- xDelete<Input*>( this -> inputs );
64- xDelete<IssmDouble>( this -> timesteps );
50+
51+ /* No need to delete vectors, they manage their own memory */
6552
6653 if (this ->current_input ) delete this ->current_input ;
6754}
@@ -70,22 +57,23 @@ TransientInput::~TransientInput(){/*{{{*/
7057/* Object virtual functions definitions:*/
7158Input* TransientInput::copy () {/* {{{*/
7259
60+ _assert_ (this ->timesteps .size ()==this ->numtimesteps );
61+ _assert_ (this ->inputs .size ()==this ->numtimesteps );
62+
7363 TransientInput* output=NULL ;
7464
7565 output = new TransientInput ();
76- output->enum_type =this ->enum_type ;
77- output->numtimesteps =this ->numtimesteps ;
78- if (this ->numtimesteps >0 ){
79- output->timesteps =xNew<IssmDouble>(this ->numtimesteps );
80- xMemCpy (output->timesteps ,this ->timesteps ,this ->numtimesteps );
81- output->inputs = xNew<Input*>(this ->numtimesteps );
82- for (int i=0 ;i<this ->numtimesteps ;i++){
83- if (this ->inputs [i]){
84- output->inputs [i] = this ->inputs [i]->copy ();
85- }
86- else {
87- output->inputs [i] = NULL ;
88- }
66+ output->enum_type = this ->enum_type ;
67+ output->numtimesteps = this ->enum_type ;
68+
69+ output->timesteps = this ->timesteps ; /* performs a Deep copy*/
70+
71+ for (int i=0 ;i<this ->numtimesteps ;i++){
72+ if (this ->inputs [i]){
73+ output->inputs .push_back (this ->inputs [i]->copy ());
74+ }
75+ else {
76+ output->inputs .push_back (NULL );
8977 }
9078 }
9179 output->parameters =this ->parameters ;
@@ -125,7 +113,11 @@ int TransientInput::Id(void){ return -1; }/*{{{*/
125113/* }}}*/
126114void TransientInput::Marshall (MarshallHandle* marshallhandle){ /* {{{*/
127115
128- bool isnull;
116+ _assert_ (this ->timesteps .size ()==this ->numtimesteps );
117+ _assert_ (this ->inputs .size ()==this ->numtimesteps );
118+
119+ bool isnull;
120+ IssmDouble time;
129121
130122 int object_enum = TransientInputEnum;
131123 marshallhandle->call (object_enum);
@@ -134,25 +126,15 @@ void TransientInput::Marshall(MarshallHandle* marshallhandle){ /*{{{*/
134126 marshallhandle->call (this ->numberofvertices_local );
135127 marshallhandle->call (this ->enum_type );
136128 marshallhandle->call (this ->numtimesteps );
137- marshallhandle->call (this ->timesteps ,numtimesteps);
138-
139- /* Allocate memory if need be*/
140- if (marshallhandle->OperationNumber () == MARSHALLING_LOAD){
141- int N = this ->numtimesteps ; _assert_ (N>=0 && N<1e6 );
142- if (N){
143- this ->inputs = xNew<Input*>(N);
144- for (int i=0 ;i<N;i++) this ->inputs [i] = NULL ;
145- }
146- else {
147- this ->inputs = NULL ;
148- }
149- }
150129
151130 /* Marshall!*/
152131 if (marshallhandle->OperationNumber ()!=MARSHALLING_LOAD){
153132 for (int i=0 ;i<this ->numtimesteps ;i++){
154133
155- // _assert_(this->inputs[i]);
134+ /* Write time for slice i*/
135+ marshallhandle->call (this ->timesteps [i]);
136+
137+ /* Write input for slice i*/
156138 isnull = false ;
157139 if (!this ->inputs [i]) isnull = true ;
158140 marshallhandle->call (isnull);
@@ -166,27 +148,39 @@ void TransientInput::Marshall(MarshallHandle* marshallhandle){ /*{{{*/
166148 }
167149 else {
168150 for (int i=0 ;i<this ->numtimesteps ;i++){
151+
152+ /* Load time for slice i*/
153+ marshallhandle->call (time);
154+ this ->timesteps .push_back (time);
155+
156+ /* Load input for slice i*/
169157 marshallhandle->call (isnull);
170158 if (!isnull){
171159 marshallhandle->call (object_enum);
172160
173161 if (object_enum==TriaInputEnum){
174162 TriaInput* triainput2=new TriaInput ();
175163 triainput2->Marshall (marshallhandle);
176- this ->inputs [i]= triainput2;
164+ this ->inputs . push_back ( triainput2) ;
177165 }
178166 else if (object_enum==PentaInputEnum){
179167 PentaInput* pentainput2=new PentaInput ();
180168 pentainput2->Marshall (marshallhandle);
181- this ->inputs [i]= pentainput2;
169+ this ->inputs . push_back ( pentainput2) ;
182170 }
183171 else {
184172 _error_ (" input " <<EnumToStringx (object_enum)<<" not supported" );
185173 }
186174 }
175+ else {
176+ this ->inputs .push_back (NULL );
177+ }
187178 }
188179 }
189180
181+ _assert_ (this ->timesteps .size ()==this ->numtimesteps );
182+ _assert_ (this ->inputs .size ()==this ->numtimesteps );
183+
190184}
191185/* }}}*/
192186int TransientInput::ObjectEnum (void ){/* {{{*/
@@ -208,35 +202,14 @@ void TransientInput::AddTriaTimeInput(IssmDouble time,int numindices,int* indice
208202 }
209203
210204 /* This is a new time step! we need to add it to the list*/
211- if (this ->numtimesteps >0 && time<this ->timesteps [this ->numtimesteps -1 ]) _error_ (" timestep values must increase sequentially, here " << this ->timesteps [this ->numtimesteps -1 ] <<" is the last step but smaller than the preceding " << time<<" \n " );
212-
213- IssmDouble *old_timesteps = NULL ;
214- Input **old_inputs = NULL ;
215- if (this ->numtimesteps > 0 ){
216- old_timesteps=xNew<IssmDouble>(this ->numtimesteps );
217- xMemCpy (old_timesteps,this ->timesteps ,this ->numtimesteps );
218- xDelete<IssmDouble>(this ->timesteps );
219- old_inputs=xNew<Input*>(this ->numtimesteps );
220- xMemCpy (old_inputs,this ->inputs ,this ->numtimesteps );
221- xDelete<Input*>(this ->inputs );
205+ if (this ->numtimesteps >0 && time<this ->timesteps [this ->numtimesteps -1 ]){
206+ _error_ (" timestep values must increase sequentially, here " << this ->timesteps [this ->numtimesteps -1 ] <<" is the last step but smaller than the preceding " << time<<" \n " );
222207 }
223208
224- this ->numtimesteps =this ->numtimesteps +1 ;
225- this ->timesteps =xNew<IssmDouble>(this ->numtimesteps );
226- this ->inputs = xNew<Input*>(this ->numtimesteps );
227-
228- if (this ->numtimesteps > 1 ){
229- xMemCpy (this ->inputs ,old_inputs,this ->numtimesteps -1 );
230- xMemCpy (this ->timesteps ,old_timesteps,this ->numtimesteps -1 );
231- xDelete (old_timesteps);
232- xDelete<Input*>(old_inputs);
233- }
234-
235- /* go ahead and plug: */
236- this ->timesteps [this ->numtimesteps -1 ] = time;
237- this ->inputs [this ->numtimesteps -1 ] = NULL ;
209+ this ->numtimesteps ++;
210+ this ->timesteps .push_back (time);
211+ this ->inputs .push_back (NULL );
238212 this ->AddTriaTimeInput (this ->numtimesteps -1 ,numindices,indices,values_in,interp_in);
239-
240213}
241214/* }}}*/
242215void TransientInput::AddPentaTimeInput (IssmDouble time,int numindices,int * indices,IssmDouble* values_in,int interp_in){/* {{{*/
@@ -252,33 +225,10 @@ void TransientInput::AddPentaTimeInput(IssmDouble time,int numindices,int* indic
252225 /* This is a new time step! we need to add it to the list*/
253226 if (this ->numtimesteps >0 && time<this ->timesteps [this ->numtimesteps -1 ]) _error_ (" timestep values must increase sequentially" );
254227
255- IssmDouble *old_timesteps = NULL ;
256- Input **old_inputs = NULL ;
257- if (this ->numtimesteps > 0 ){
258- old_timesteps=xNew<IssmDouble>(this ->numtimesteps );
259- xMemCpy (old_timesteps,this ->timesteps ,this ->numtimesteps );
260- xDelete<IssmDouble>(this ->timesteps );
261- old_inputs=xNew<Input*>(this ->numtimesteps );
262- xMemCpy (old_inputs,this ->inputs ,this ->numtimesteps );
263- xDelete<Input*>(this ->inputs );
264- }
265-
266- this ->numtimesteps =this ->numtimesteps +1 ;
267- this ->timesteps =xNew<IssmDouble>(this ->numtimesteps );
268- this ->inputs = xNew<Input*>(this ->numtimesteps );
269-
270- if (this ->numtimesteps > 1 ){
271- xMemCpy (this ->inputs ,old_inputs,this ->numtimesteps -1 );
272- xMemCpy (this ->timesteps ,old_timesteps,this ->numtimesteps -1 );
273- xDelete (old_timesteps);
274- xDelete<Input*>(old_inputs);
275- }
276-
277- /* go ahead and plug: */
278- this ->timesteps [this ->numtimesteps -1 ] = time;
279- this ->inputs [this ->numtimesteps -1 ] = NULL ;
228+ this ->numtimesteps ++;
229+ this ->timesteps .push_back (time);
230+ this ->inputs .push_back (NULL );
280231 this ->AddPentaTimeInput (this ->numtimesteps -1 ,numindices,indices,values_in,interp_in);
281-
282232}
283233/* }}}*/
284234void TransientInput::AddTriaTimeInput (int step,int numindices,int * indices,IssmDouble* values_in,int interp_in){/* {{{*/
@@ -321,7 +271,7 @@ void TransientInput::GetAllTimes(IssmDouble** ptimesteps,int* pnumtimesteps){/*{
321271
322272 if (ptimesteps){
323273 *ptimesteps=xNew<IssmDouble>(this ->numtimesteps );
324- xMemCpy (*ptimesteps,this ->timesteps ,this ->numtimesteps );
274+ xMemCpy (*ptimesteps,this ->timesteps . data () ,this ->numtimesteps );
325275 }
326276 if (pnumtimesteps){
327277 *pnumtimesteps = this ->numtimesteps ;
@@ -478,10 +428,10 @@ void TransientInput::SetCurrentTimeInput(IssmDouble time){/*{{{*/
478428
479429 /* Figure step out*/
480430 int offset, prevoffset;
481- if (!binary_search (&offset,time,this ->timesteps , this ->numtimesteps )){
431+ if (!binary_search (&offset, time, this ->timesteps . data (), this ->numtimesteps )){
482432 _error_ (" Input not found (is TransientInput sorted ?)" );
483433 }
484- if (!binary_search (&prevoffset,reCast<IssmDouble>(time-dt),this ->timesteps , this ->numtimesteps )){
434+ if (!binary_search (&prevoffset, reCast<IssmDouble>(time-dt), this ->timesteps . data (), this ->numtimesteps )){
485435 _error_ (" Input not found (is TransientInput sorted ?)" );
486436 }
487437
@@ -564,9 +514,9 @@ void TransientInput::SetAverageAsCurrentTimeInput(IssmDouble start_time,IssmDoub
564514 int found,start_offset,end_offset,input_offset;
565515
566516 /* go through the timesteps, and grab offset for start and end*/
567- found= binary_search (&start_offset,start_time,this ->timesteps , this ->numtimesteps );
517+ found = binary_search (&start_offset, start_time,this ->timesteps . data (), this ->numtimesteps );
568518 if (!found) _error_ (" Input not found (is TransientInput sorted ?)" );
569- found= binary_search (&end_offset,end_time,this ->timesteps , this ->numtimesteps );
519+ found = binary_search (&end_offset, end_time, this ->timesteps . data (), this ->numtimesteps );
570520 if (!found) _error_ (" Input not found (is TransientInput sorted ?)" );
571521
572522 if (start_offset==-1 ){
@@ -640,6 +590,7 @@ void TransientInput::SetAverageAsCurrentTimeInput(IssmDouble start_time,IssmDoub
640590 }
641591 _assert_ (dtsum>0 );
642592 durinv=1 ./dtsum;
593+
643594 /* Integration done, now normalize*/
644595 switch (averaging_method){
645596 case 0 : // Arithmetic mean
@@ -656,18 +607,20 @@ void TransientInput::SetAverageAsCurrentTimeInput(IssmDouble start_time,IssmDoub
656607 }
657608}/* }}}*/
658609IssmDouble TransientInput::GetTimeByOffset (int offset){/* {{{*/
610+
659611 if (offset<0 ) offset=0 ;
660612 _assert_ (offset<this ->numtimesteps );
661613 return this ->timesteps [offset];
614+
662615}
663616/* }}}*/
664617int TransientInput::GetTimeInputOffset (IssmDouble time){/* {{{*/
665618
666- int offset;
667-
668619 /* go through the timesteps, and figure out which interval we
669- * *fall within. Then interpolate the values on this interval: */
670- int found=binary_search (&offset,time,this ->timesteps ,this ->numtimesteps );
620+ * fall within. Then interpolate the values on this interval: */
621+
622+ int offset;
623+ int found = binary_search (&offset, time, this ->timesteps .data (), this ->numtimesteps );
671624 if (!found) _error_ (" Input not found (is TransientInput sorted ?)" );
672625
673626 return offset;
0 commit comments