source: palm/trunk/UTIL/chemistry/gasphase_preproc/kpp4palm/src/create_kpp_module.C @ 3797

Last change on this file since 3797 was 3797, checked in by forkel, 3 years ago

Modifications for OpenMP version by Klaus Ketelsen

File size: 25.4 KB
Line 
1
2// ############################################################################
3//
4//     create_kpp_module
5//
6//     create scalar code from .f90 sources created by KPP
7//
8//     COPYRIGHT Klaus Ketelsen and MPI-CH   April 2007
9//
10// ############################################################################
11//
12//Current revisions:
13//------------------
14//
15//
16//Former revisions:
17//-----------------
18//$Id: create_kpp_module.C 3453 2018-10-30 13:21:51Z forkel $
19// OpenMP version    (15.03.2019, ketelsen)
20//
21// Added vector switch Kacc,Krej,IERRV, Commented add_line for istatf,    (05.03.2019, forkel)
22//      added ,pe after ierr_u,         
23//
24// Added create_set_cs and cs_mech and get_mechanismname in module_header (05.03.2019, forkel)
25//
26// exclude kco_compress from handling by global_variables2vector (30.10.2018, forkel)
27//
28// Added  automatic line with mechanism name (read mech_list) (25.09.2018, forkel)
29//
30// Added  vl_glo = size(tempi,1) (20.09.2018, forkel)
31//
32// Removed creation of fill_ Subroutine and creation of calls thereof (18.09.2018, ketelsen)
33//
34// Fix in order not to loose the values of qvap and fakt (12.09.2018, forkel)
35//
36// Bug fixes: moved kppi.add_line("    CALL initialize after fakt = fakti(is)
37// Deleted definition of qvap,fakt in create_kpp_integrate again (03.09.2018, forkel)
38//
39// Changes for vector mode (edit_WAXPY, edit_FunTemplate, edit_JacTemplate,
40// some cleanup of comments, various changes in create_kpp_integrate) (July 2018, ketelsen)
41//
42//
43// Added qvap and fakt                            (June 2018, forkel)
44// --> Change in module_header: qvap, fakt added  (June 2018, forkel)
45//
46// re-established original uppercase/lowercase     (June 2018, forkel)
47// --> Change in module_header: reset case in  Initialize, Integrate, and
48//     Update_rconst                               (June 2018, forkel)
49//
50// Removed preprocessor directive __chem again (2017-09-14, forkel)
51//
52// Added phot                                                 (2017-09-14, forkel)
53// --> Change in module_header: Variables for photolyis added (2017-09-14, forkel)
54//
55// change of some output to lowercase with uppercase Fortran (2017, forkel)
56//
57// Intial version of KP4 adapted to PALM                      (Nov. 2016, ketelsen)
58//
59#include <stdio.h>
60// stdlib is necessary to define getenv:
61#include <stdlib.h>
62
63#include "create_kpp_module.h"
64#include "utils.h"
65
66void create_kpp_module::do_work (string s) {
67   vector<fortran_file>::iterator  it;
68   vector<string>::iterator        ic;
69   vector<Vvar>::iterator          iv;
70
71   expand_decomp                   exp_de;
72
73   prefix = s;
74   module_name = prefix;
75
76   cout << "Create " << module_name << " from kpp Fortran sources" <<endl;
77   cout << "Vector mode " << kpp_switches.is_vector() <<endl;
78   cout << "De_indexing " << kpp_switches.de_indexing() <<endl;
79
80   create_fortran_files_and_read();
81
82// Generate first module lines
83
84     string first_line="MODULE " + module_name;
85   mz_kpp.add_line(first_line);
86   mz_kpp.add_line(" ");
87
88//    string e5_line = first_line +"_e5";
89//    e5_kpp.add_line(e5_line);
90//    e5_line = "  USE             " + module_name;
91//    e5_kpp.add_line(e5_line);
92//    e5_kpp.add_line(" ");
93
94// edit include files
95
96   for(it=kpp_includes.begin();it!=kpp_includes.end();it++) {
97     it->edit_inc(header_variables);
98
99//   Create variable Species list and vector variable list
100
101     if(it->get_name() == module_name + "_Parameters") {
102       it->create_species_list(species_list);
103     }
104     if(it->get_name() == module_name + "_Global") {
105       it->vector_variable_list(Vvar_list);
106     }
107   }
108
109// Prepare expansion of decomposition subroutine
110
111   if(kpp_switches.de_indexing () > 0 ) {
112     exp_de.create_sparse_info (kpp_includes, module_name);
113   }
114
115// edit FORTRAN files
116
117   for(it=kpp_files.begin();it!=kpp_files.end();it++) {
118     it->edit_fortran ();
119   }
120
121// Generate a list of single subroutines from kpp-files
122// kpp files are modules containing several subroutines
123
124   copy_files_to_subroutines ();
125
126// All header_variables to include list
127
128   kpp_includes.push_back(header_variables);
129
130// Create decomposition subroutine
131   if(kpp_switches.de_indexing () > 0 ) {
132     exp_de.create_routine (kpp_subroutines);
133   }
134
135   if(kpp_switches.is_vector()) {
136
137   cout << "##### Hier kpp_switches.is_vector          " <<endl;
138//   Change header section
139     for(it=kpp_includes.begin();it!=kpp_includes.end();it++) {
140       it->edit_inc_vec(global_variable_list);
141     }
142
143//   Change global variables to vector (except for kp4_compress, which has already the right form)
144     
145     for(it=kpp_subroutines.begin();it!=kpp_subroutines.end();it++) {
146       if(it->get_name() != "kco_compress" ) {
147         it->global_variables2vector (global_variable_list);
148       }
149     }
150
151//   Edit individual subroutines
152
153     for(it=kpp_subroutines.begin();it!=kpp_subroutines.end();it++) {
154       if(it->get_name() == "KppDecomp") {
155         it->edit_KppDecomp();
156       }
157       if(it->get_name() == "KppSolve") {
158         it->edit_KppSolve();
159       }
160       if(it->get_name() == "Jac_SP" ) {
161         it->edit_Jac_SP();
162       }
163       if(it->get_name() == "Fun" ) {
164         it->edit_Fun();
165       }
166       if(it->get_name() == "WAXPY" ) {
167         it->edit_WAXPY();
168       }
169       if(it->get_name() == "FunTemplate" ) {
170         it->edit_FunTemplate();
171       }
172       if(it->get_name() == "JacTemplate" ) {
173         it->edit_JacTemplate();
174       }
175     }
176   }
177
178// Update_RCONST has to be changed also in scalar mode
179
180   for(it=kpp_subroutines.begin();it!=kpp_subroutines.end();it++) {
181     if(it->get_name() == "Update_RCONST") {
182       it->edit_Update_RCONST(Vvar_list);
183     }
184
185     if(it->get_name() == "Initialize") {
186       it->edit_Initialize(Vvar_list);
187     }
188
189   }
190
191// Add Solver template to subroutine list
192   if(kpp_switches.is_vector()) {
193     add_solver_to_subroutine_list ();
194   }
195
196// The module header will be taken from ../templates/module_header.
197// Please edit if header has to be changed.
198
199   generate_module_header();
200
201// create_set_cs
202   create_set_cs();
203
204// Create kpp_integrate subroutine (chem_gasphase_integrate) for skalar and vector mode
205
206   create_kpp_integrate();
207// Copy include files
208
209   for(it=kpp_includes.begin();it!=kpp_includes.end();it++) {
210     it->copy_to_MZ_KPP(mz_kpp);
211   }
212
213   mz_kpp.add_line(" ");
214   mz_kpp.add_line("! Interface Block ");
215   mz_kpp.add_line(" ");
216   for(it=kpp_subroutines.begin();it!=kpp_subroutines.end();it++) {
217     string          buf;
218
219     string prefix = "  ";
220     for(ic=interface_ignore.begin();ic!=interface_ignore.end();ic++) {
221       if(it->get_name() == *ic) {
222         prefix = "!interface not working  ";
223         break;
224       }
225     }
226
227     buf = prefix + "interface            " + it->get_name() ;
228     mz_kpp.add_line(buf);
229     buf = prefix + "  module procedure   " + it->get_name();
230     mz_kpp.add_line(buf);
231     buf = prefix + "end interface        " + it->get_name();
232     mz_kpp.add_line(buf);
233     mz_kpp.add_line(" ");
234   }
235
236   mz_kpp.add_line(" ");
237
238// Declare variables THREADPRIVATE for OpenMP version
239
240   mz_kpp.add_line("  ! OpenMP directives generated by kp4 ");
241   mz_kpp.add_line(" ");
242   mz_kpp.add_line("  !$OMP THREADPRIVATE (vl,vl_glo,is,ie,data_loaded)");
243   mz_kpp.add_line("  !$OMP THREADPRIVATE (c,var,fix,rconst,time,temp,stepmin,cfactor)");
244   mz_kpp.add_line("  !$OMP THREADPRIVATE (qvap,fakt,cs_mech,a,icntrl,rcntrl)");
245   mz_kpp.add_line(" ");
246   if(kpp_switches.is_vector()) {
247      mz_kpp.add_line("  ! Vector mode Only ");
248          mz_kpp.add_line(" ");
249          mz_kpp.add_line("  !$OMP THREADPRIVATE (kacc,krej,ierrv)");
250          mz_kpp.add_line("  !$OMP THREADPRIVATE (kpoints,kpoints_SAVE,index_org,done_check,index_step,cell_done)");
251          mz_kpp.add_line("  !$OMP THREADPRIVATE (f_done,kacc_done,krej_done,ierr_done,compress_done)");
252          mz_kpp.add_line(" ");
253   }
254
255
256// Copy FORTRAN subroutines to mz_kpp
257
258   mz_kpp.add_line(" CONTAINS");
259   
260   for(it=kpp_subroutines.begin();it!=kpp_subroutines.end();it++) {
261     mz_kpp.add_line(" ");
262     it->copy_to_MZ_KPP(mz_kpp);
263   }
264
265// Finish module
266
267   string last_line="end module " + module_name;
268   mz_kpp.add_line("");
269   mz_kpp.add_line(last_line);
270
271// Write the complete module to file: mz_kpp.f
272
273   write_module_file();
274
275   return;
276}
277
278void create_kpp_module::create_fortran_files_and_read() {
279
280   string                          name;
281   ifstream                        in,in_c,in_b,in_i;
282   fortran_file                    f_file;
283   vector<fortran_file>::iterator  it;
284
285// Open file with list of FORTRAN routines
286
287   in.open("file_list");
288   if( !in ) {
289      cout << "cannot open " << endl; my_abort("file_list");
290   }
291   
292// Create kpp_fortran routines
293   while ( 1 ) {
294     in >> name;
295     if( in.eof() ) break;
296     if( in.bad() ) my_abort("ERROR_READ_1");
297     f_file.set_name(name);
298     kpp_files.push_back(f_file);
299   }
300   in.close();
301
302// Read FORTRAN code
303
304   for(it=kpp_files.begin();it!=kpp_files.end();it++) {
305     it->read();
306   }
307
308// Open file with list of include files
309
310   in_c.open("include_list");
311   if( !in_c ) {
312      cout << "cannot open " << endl; my_abort("include_list");
313   }
314
315// Create kpp_includes vector
316   while ( 1 ) {
317     in_c >> name;
318     if( in_c.eof() ) break;
319     if( in_c.bad() ) my_abort("ERROR_READ_3");
320     f_file.set_name(name);
321     kpp_includes.push_back(f_file);
322   }
323   in_c.close();
324
325// Read include files
326
327   for(it=kpp_includes.begin();it!=kpp_includes.end();it++) {
328     it->read();
329   }
330
331// Read Ignore list
332
333   in_i.open("interface_ignore_list");
334   if( !in_i ) {
335      cout << "cannot open " << endl; my_abort("include_list");
336   }
337
338// Create kpp_includes vector
339   while ( 1 ) {
340     in_i >> name;
341     if( in_i.eof() ) break;
342     if( in_i.bad() ) my_abort("ERROR_READ_4");
343     interface_ignore.push_back(name);
344   }
345   in_c.close();
346
347}
348
349void create_kpp_module::copy_files_to_subroutines () {
350   string                          name;
351   ifstream                        in;
352   fortran_file                    s_file;
353   vector<fortran_file>::iterator  it;
354
355// Open file with list of FORTRAN routines
356
357   in.open("subroutine_list");
358   if( !in ) {
359      cout << "cannot open " << endl; my_abort("subroutine_list");
360   }
361
362// Create vector kpp_subroutines
363
364   while ( 1 ) {
365     in >> name;
366     if( in.eof() ) break;
367     if( in.bad() ) my_abort("ERROR_READ_S1");
368     s_file.set_name(name);
369     kpp_subroutines.push_back(s_file);
370   }
371   in.close();
372
373   header_variables.add_line(" ");
374   header_variables.add_line("!  variable definations from  individual module headers ");
375   header_variables.add_line(" ");
376
377//  Loop over all FORTRAN Files
378
379   for(it=kpp_files.begin();it!=kpp_files.end();it++) {
380     it->copy_to_subroutine_vector(kpp_subroutines, header_variables);
381   }
382}
383
384void create_kpp_module::add_solver_to_subroutine_list () {
385   fortran_file                    s_file;
386
387   string solver_name = getenv("KPP_SOLVER");
388   cout << "KPP_SOLVER " <<solver_name <<endl;
389   
390   s_file.set_name(solver_name);
391   s_file.read();
392   kpp_subroutines.push_back(s_file);
393
394   return;
395}
396
397void create_kpp_module::generate_module_header() {
398
399   string                          buf;
400   ifstream                        in;
401   ifstream                        in_e5;
402   program_line                    line;
403   vector<fortran_file>::iterator  it;
404   char                            distr[2];
405   string                          diline;
406
407// Read mechanism from mech_list
408
409   in.open("mech_list");
410   if( !in ) {
411      cout << "cannot open " << endl; my_abort("mech_list");
412   }
413
414   while ( 1 ) {
415     getline (in, buf);
416     if( in.eof() ) break;
417     if( in.bad() ) my_abort("ERROR_READ_4");
418     line.set_line(buf);
419     mz_kpp.add_line(line);
420   }
421   in.close();
422
423
424// Read Modul Header from file $MZ_KPP_HOME/templates/module_header
425
426   in.open("module_header");
427   if( !in ) {
428      cout << "cannot open " << endl; my_abort("module_header");
429   }
430
431   while ( 1 ) {
432     getline (in, buf);
433     if( in.eof() ) break;
434     if( in.bad() ) my_abort("ERROR_READ_4");
435     line.set_line(buf);
436     mz_kpp.add_line(line); 
437   }
438   mz_kpp.add_line("                                                                 "); 
439   mz_kpp.add_line("! Variables used for vector mode                                 "); 
440   mz_kpp.add_line("                                                                 "); 
441   if(kpp_switches.is_vector()) {
442       mz_kpp.add_line("  logical,parameter          :: L_VECTOR = .TRUE.             ");
443   } else {
444       mz_kpp.add_line("  logical,parameter          :: L_VECTOR = .FALSE.            ");
445   }
446//  mz_pj_20070531+
447   sprintf(distr,"%i",kpp_switches.de_indexing());
448   diline = distr ;
449   mz_kpp.add_line("  integer,parameter          :: I_LU_DI = " + diline );
450//  mz_pj_20070531-
451
452   mz_kpp.add_line("  integer,parameter          :: VL_DIM = " 
453                 + kpp_switches.get_vector_length() ); 
454   mz_kpp.add_line("  integer                     :: vl                              "); 
455   mz_kpp.add_line("                                                                 "); 
456   mz_kpp.add_line("  integer                     :: VL_glo                          "); 
457   mz_kpp.add_line("  integer                     :: is,ie                           "); 
458   mz_kpp.add_line("                                                                 "); 
459   mz_kpp.add_line("                                                                 "); 
460   if(kpp_switches.is_vector()) {
461      mz_kpp.add_line("  integer, dimension(VL_dim)   :: Kacc,Krej                       "); 
462      mz_kpp.add_line("  integer, dimension(VL_dim)   :: IERRV                           "); 
463   }
464   mz_kpp.add_line("  logical                     :: data_loaded = .false.             "); 
465   if(kpp_switches.is_vector()) {
466          mz_kpp.add_line("  REAL(dp),POINTER,DIMENSION(:,:),CONTIGUOUS    :: var           ");
467   } else {
468      mz_kpp.add_line("  REAL(dp),POINTER,DIMENSION(:),CONTIGUOUS    :: var             ");
469   }
470   in.close();
471
472   return;
473}
474
475void create_kpp_module::write_module_file() {
476   ofstream                    out;
477   ofstream                    out_e5;
478
479   string out_file  = "kk_kpp.f90";
480   out.open(out_file.c_str(), ios::out);
481   if( !out ) {
482      cout << "cannot open " << endl; my_abort(out_file);
483   }
484
485   mz_kpp.write_file (out);
486
487   out.close();
488   
489
490   return;
491}
492
493void create_kpp_module::create_set_cs() {
494   fortran_file          kppi;         
495   vector<Vvar>::iterator               iv;
496   string                               xline;
497     
498   string                          buf;
499   ifstream                        in;
500   program_line                    line;
501
502   kppi.set_name("get_mechanismname");
503   kppi.add_line("SUBROUTINE get_mechanismname                                        ");
504   kppi.add_line("                                                                    ");
505   kppi.add_line("  IMPLICIT NONE                                                     ");
506// Read mechanism from set_cm
507// Tis got an own own subroutine to aviod being called at each timestep
508
509   in.open("set_cm");
510   if( !in ) {
511      cout << "cannot open " << endl; my_abort("set_cm");
512   }
513
514   while ( 1 ) {
515     getline (in, buf);
516     if( in.eof() ) break;
517     if( in.bad() ) my_abort("ERROR_READ_4");
518     line.set_line(buf);
519     kppi.add_line(line);
520   }
521   in.close();
522
523   kppi.add_line("                                                                    ");
524   kppi.add_line("  return                                                            ");
525   kppi.add_line("END SUBROUTINE get_mechanismname                                    ");
526   kppi.add_line("                                                                    ");
527   kpp_subroutines.push_back(kppi);
528
529   return;
530}
531
532
533void create_kpp_module::create_kpp_integrate() {
534   fortran_file          kppi;
535   vector<Vvar>::iterator               iv;
536   string                               xline;
537
538
539   kppi.set_name("chem_gasphase_integrate");
540
541   kppi.add_line("SUBROUTINE chem_gasphase_integrate (time_step_len, conc, tempi, qvapi, fakti, photo, ierrf, xnacc, xnrej, istatus, l_debug, pe, icntrl_i, rcntrl_i )  ");
542   kppi.add_line("                                                                    ");
543   kppi.add_line("  IMPLICIT NONE                                                     ");
544   kppi.add_line("                                                                    ");
545
546   kppi.add_line("  REAL(dp), INTENT(IN)                   :: time_step_len           ");
547   kppi.add_line("  REAL(dp),  DIMENSION(:,:),  INTENT(INOUT) :: conc                    ");
548   kppi.add_line("  REAL(dp),  DIMENSION(:,:),  INTENT(IN)    :: photo                   ");
549   kppi.add_line("  REAL(dp),  DIMENSION(:),  INTENT(IN)      :: tempi                   ");
550   kppi.add_line("  REAL(dp),  DIMENSION(:),  INTENT(IN)      :: qvapi                   ");
551   kppi.add_line("  REAL(dp),  DIMENSION(:),  INTENT(IN)      :: fakti                   ");
552   kppi.add_line("  INTEGER,  INTENT(OUT), OPTIONAL        :: ierrf(:)                ");
553   kppi.add_line("  INTEGER,  INTENT(OUT), OPTIONAL        :: xNacc(:)                ");
554   kppi.add_line("  INTEGER,  INTENT(OUT), OPTIONAL        :: xNrej(:)                ");
555   kppi.add_line("  INTEGER,  INTENT(INOUT), OPTIONAL      :: istatus(:)              ");
556   kppi.add_line("  INTEGER,  INTENT(IN), OPTIONAL         :: PE                      ");
557   kppi.add_line("  LOGICAL,  INTENT(IN), OPTIONAL         :: l_debug                 ");
558   kppi.add_line("  INTEGER,  DIMENSION(nkppctrl),INTENT(IN), OPTIONAL  :: icntrl_i         ");
559   kppi.add_line("  REAL(dp), DIMENSION(nkppctrl),INTENT(IN), OPTIONAL  :: rcntrl_i         ");
560   kppi.add_line("                                                                    ");
561   kppi.add_line("  INTEGER                                 :: k   ! loop variable     ");
562   kppi.add_line("  REAL(dp)                                :: dt                      ");
563   kppi.add_line("  integer, dimension(20)                 :: istatus_u               ");
564   kppi.add_line("  integer                                :: ierr_u                  ");
565// kppi.add_line("  integer                                :: istatf                  ");
566   kppi.add_line("  integer                                :: vl_dim_lo               ");
567   kppi.add_line("                                                                    ");
568   kppi.add_line("                                                                    ");
569   kppi.add_line("  if (present (istatus) )   istatus = 0                             ");
570   kppi.add_line("  if (present (icntrl_i) )  icntrl  = icntrl_i                      ");
571   kppi.add_line("  if (present (rcntrl_i) )  rcntrl  = rcntrl_i                      ");
572   kppi.add_line("                                                                    ");
573   if(kpp_switches.is_vector()) {
574      kppi.add_line("  var => c(:,1:nvar)                                                  ");
575   } else {
576      kppi.add_line("  var => c(1:nvar)                                                  ");
577   }
578   kppi.add_line("                                                                    ");
579   kppi.add_line("  vl_glo = size(tempi,1)                                            ");
580   kppi.add_line("                                                                    ");
581   kppi.add_line("  vl_dim_lo = VL_DIM                                                ");
582   kppi.add_line("  DO k=1,VL_glo,vl_dim_lo                                           ");
583   kppi.add_line("    is = k                                                          ");
584   kppi.add_line("    ie = min(k+vl_dim_lo-1,VL_glo)                                  ");
585   kppi.add_line("    vl = ie-is+1                                                    ");
586
587   kppi.add_line("                                                                    ");
588   if(kpp_switches.is_vector()) {
589     kppi.add_line("    C(1:vl,:) = Conc(is:ie,:)                                     ");
590   } else {
591     kppi.add_line("    C(:) = Conc(is,:)                                             ");
592   }
593
594   kppi.add_line("                                                                    ");
595   if(kpp_switches.is_vector()) {
596     kppi.add_line("    temp(1:vl) = tempi(is:ie)                                     ");
597   } else {
598     kppi.add_line("    temp = tempi(is)                                              ");
599   }
600   kppi.add_line("                                                                    ");
601   if(kpp_switches.is_vector()) {
602     kppi.add_line("    qvap(1:vl) = qvapi(is:ie)                                     ");
603   } else {
604     kppi.add_line("    qvap = qvapi(is)                                              ");
605   }
606   kppi.add_line("                                                                    ");
607   if(kpp_switches.is_vector()) {
608     kppi.add_line("    fakt(1:vl) = fakti(is:ie)                                     ");
609   } else {
610     kppi.add_line("    fakt = fakti(is)                                              ");
611   }
612
613   kppi.add_line("                                                                    ");
614   kppi.add_line("    CALL initialize                                                 ");
615
616   kppi.add_line("                                                                    ");
617   if(kpp_switches.is_vector()) {
618     kppi.add_line("    phot(1:vl,:) = photo(is:ie,:)                                     ");
619   } else {
620     kppi.add_line("    phot(:) = photo(is,:)                                             ");
621   }
622   kppi.add_line("                                                                    ");
623   kppi.add_line("    CALL update_rconst                                              ");
624   kppi.add_line("                                                                    ");
625   kppi.add_line("    dt = time_step_len                                              ");
626   kppi.add_line("                                                                    ");
627   kppi.add_line("    ! integrate from t=0 to t=dt                                    ");
628   kppi.add_line("    CALL integrate(0._dp, dt, icntrl, rcntrl, istatus_u = istatus_u, ierr_u=ierr_u)");
629   kppi.add_line("                                                                    ");
630   kppi.add_line("                                                                    ");
631   if(kpp_switches.is_vector()) {
632     kppi.add_line("    Conc(is:ie,:) = C(1:VL,:)                                     ");
633   } else {
634     kppi.add_line("   IF (PRESENT(l_debug) .AND. PRESENT(PE)) THEN                       ");
635     kppi.add_line("      IF (l_debug) CALL error_output(Conc(is,:),ierr_u,pe)            ");
636     kppi.add_line("   ENDIF                                                              ");
637     kppi.add_line("                                                                      ");
638     kppi.add_line("    Conc(is,:) = C(:)                                                 ");
639   }
640
641   kppi.add_line("                                                                    ");
642   kppi.add_line("    ! Return Diagnostic Information                                 ");
643   kppi.add_line("                                                                    ");
644   if(kpp_switches.is_vector()) {
645     kppi.add_line("    if(Present(ierrf))    ierrf(is:ie) = IERRV(1:VL)              ");
646     kppi.add_line("    if(Present(xNacc))    xNacc(is:ie) = Kacc(1:VL)               ");
647     kppi.add_line("    if(Present(xNrej))    xNrej(is:ie) = Krej(1:VL)               ");
648   } else {
649     kppi.add_line("    if(Present(ierrf))    ierrf(is) = IERR_U                      ");
650     kppi.add_line("    if(Present(xNacc))    xNacc(is) = istatus_u(4)                ");
651     kppi.add_line("    if(Present(xNrej))    xNrej(is) = istatus_u(5)                ");
652   }
653   kppi.add_line("                                                                    ");
654   kppi.add_line("    if (present (istatus) )  then                                   ");
655   if(kpp_switches.is_vector()) {
656     kppi.add_line("      istatus(4) =   istatus(4) + sum(Kacc(1:VL))                  ");
657     kppi.add_line("      istatus(5) =   istatus(5) + sum(Krej(1:VL))                  ");
658     kppi.add_line("      istatus(3) =   istatus(4) + istatus(5)                       ");
659     kppi.add_line("      istatus(6) =   istatus(6) + istatus_u(6)                     ");
660     kppi.add_line("      istatus(7) =   istatus(7) + istatus_u(7)                     ");
661   } else {
662     kppi.add_line("      istatus(1:8) = istatus(1:8) + istatus_u(1:8)                 ");
663   }
664   kppi.add_line("    end if                                                          ");
665   kppi.add_line("                                                                    ");
666   kppi.add_line("  END DO                                                            ");
667   kppi.add_line(" ");
668
669   kppi.add_line("                                                                    ");
670   kppi.add_line("! Deallocate input arrays                                           ");
671   kppi.add_line("                                                                    ");
672   for(iv=Vvar_list.begin();iv!=Vvar_list.end();iv++) {
673//     kppi.add_line("  if (allocated("+ iv->name +"))   deallocate("+ iv->name +" )    ");
674   }
675
676   kppi.add_line("                                                                    ");
677   kppi.add_line("  data_loaded = .false.                                             ");
678   kppi.add_line("                                                                    ");
679   kppi.add_line("  return                                                            ");
680   kppi.add_line("END SUBROUTINE chem_gasphase_integrate                              ");
681
682//   e5_subroutines.push_back(kppi);
683   kpp_subroutines.push_back(kppi);
684
685   return;
686}
687
Note: See TracBrowser for help on using the repository browser.