Creating a DiscreteFunction from values stored in a file

• Motivational Examples
  • Generate a log-normal distribution using external sequential Gaussian simulation software and save it to a data file. Read this distribution into Sundance and use it to simulate groundwater flow through a porous medium with spatially-varying hydraulic conductivity.
  • Generate an initial condition (IC) for a Sundance simulation using external software and save it to a data file. Read the IC into Sundance and use it to build an expression associated with the mesh that can be passed to a solver.

• Assumptions made in ReadField.cpp

  • Values of a field variable at each node in the mesh are stored in a text file with some format style. The example ReadField.cpp assumes a format whereby line 1 of the text file lists the total number of data values to be read in (N) and lines 2 through N+1 contain the data with one value per line. The user can modify the section of code

    /* fieldData.dat is the file containing the field values and nNodes is the number of values in the file */
       ifstream is("../../examples-tutorial/fieldData.dat");
       int nNodes;
       is >> nNodes;
    
       /* check for consistency between the mesh and the data */
       TEUCHOS_TEST_FOR_EXCEPTION(mesh.numCells(0) != nNodes, RuntimeError,
           "number of nodes in data file fieldData.dat is " 
           << nNodes << " but number of nodes in mesh is "
           << mesh.numCells(0));
          
          /* read the data, putting each entry into its correct place in
     the vector as indexed by the dof number, NOT the node number.
     Here I'm assuming one field value per node. */
    
       Array<int> dofs(1);
       double fVal;
       for (int i=0; i<nNodes; i++)
      {
        /* look up the dof for the 0-th function on this node */
        dofMap->getDOFsForCell(0, i, 0, dofs);
        int dof = dofs[0];
        is >> fVal;
        vec.setElement(dof, fVal);
      }
    

    which reads the data from the file to account for any variations in file formatting.

  • The field values are ordered so that they correspond to how the mesh nodes are ordered.

  • The field variable is a scalar, so the function ID in getDOFsForCell() is zero. Here is an example of how to deal with vector-valued field variables, assuming that the function values have already been saved into

    RefCountPtr<Array<Array<double> > > fVal;
    

    that has dimension (number of Nodes in mesh) x (number of elements in field variable vector). The following code fragment assumes that the field variable is a 2D vector, but the extension to higher dimensions is similar.

     for (int i = 0; i < (*fVal).size(); i++)
        {
          map->getDOFsForCell(0, i, 0, dofs);
          vec.setElement(dofs[0], (*fVal)[i][0]);
          map->getDOFsForCell(0, i, 1, dofs);
          vec.setElement(dofs[0], (*fVal)[i][1]);
        }
    

  • The code is serial. Thus when calling getDOFsForCell(), the cell local ID (cellLID_) of the current function value is i because the function value is the i-th value in the list (recall that the numbering starts at zero).

  The values are read from a data file in nodal order, then the degrees of freedom (DOF) map is consulted. The field values are reordered into a DiscreteFunction associated with the DOF in the mesh rather than the nodes, as this is how Sundance stores everything under the hood. At this point the field variable can be used within a weak form and integrated over any subset of the mesh. If one needs to access individual elements of a vector-valued field in the weak form, the square-bracket notation is used (i.e., u0[0] from the ReadField.cpp example).