Pseudo-code

The pseudo-code below shows how the ChemPlugin API might look when implemented to create a reactive transport simulator. The ChemPlugin member functions are in bold face. In actual use, you would include provision for I/O, boundary conditions, heterogeneity, multidimensionality, and so on, depending on your application.

#include <ChemPlugin.h>

void main() {
   // Create instances.
   ChemPlugin *cpi = new ChemPlugin[ninst];

   // Configure and initialize instances.
   char *config = "swap Quartz for SiO2(aq); 10 mmol/kg Quartz; …";
   for (int i=0; i<ninst; i++) {
      cpi[i].Config(config);
      cpi[i].Initialize();
   }

   // Cross-link instances.
   for (int i=1; i<ninst; i++) {
       CpiLink link1 = cpi[i].Link(cpi[i-1]);
       link1.FlowRate(flow);
       link1.Transmissivity(trans);
   }

   // Time marching loop.
   while (Time < Time_end) {
   
      // Figure stable time step.
      double deltat = 1e99;
      for (int i=0; i<ninst; i++)
         deltat = min(deltat, cpi[i].ReportTimeStep());

      // Advance time step, accounting for transport and reaction.
      for (int i=0; i<ninst; i++) cpi[i].AdvanceTimeStep(deltat);
      for (int i=0; i<ninst; i++) cpi[i].AdvanceTransport();
      for (int i=0; i<ninst; i++) cpi[i].AdvanceChemical();

      // Retrieve results.
      double qtz[ninst], perm[ninst];
      for (int i=0; i<ninst; i++)
         qtz[i] = cpi[i].Report1("saturation Quartz");
      for (int i=0; <ninst; i++)
         perm[i] = cpi[i].Report1("permeability", "darcy");
   }
}