Statistics Module

Two classes, WENOStatistics [bad link?] and StatParser [bad link?], have been created to incorporate quite flexible statistics and sampling gathering. Of these classes, only the first is directly exposed to the user.

Source code specification

In order to use this class (WENOStatistics) the user must make the following modification to the Problem.h. Note that statistics gathering is not set by default in the current solver, so a little of C++ must be done.

• include

                    #include "AMRInterpolation.h"
                    #include "WENOStatistics.h"
  

• derive your own solver class by defining

                    #define OWN_AMRSOLVER
  

• include the following protected variables in SolverSpecific

                    interpolation_type* _Interpolation; 
                    stat_type* _Stats;
  

• add the following to the bottom of the SolverSpecific constructor

                    _Interpolation = new interpolation_type();
                    _Stats = new stat_type(_Interpolation, (output_type*)_FileOutput);
  

• same for destructor

                    delete _Interpolation;
                    delete _Stats; 
  

• register a parser key with register_at(...)

                    if (_Stats) _Stats->register_at(base::LocCtrl, ""); 
  

• add the following to SetupData()

                    _Interpolation->SetupData(base::PGH(), base::NGhosts());
                    _Stats->SetupData(base::PGH(), base::NGhosts(), base::shape, base::geom);
  

• add an Evaluate call to Advance

                    _Stats->Evaluate(base::t, base::U(), base::Work());
  

• add a hook to Initialize_

                    _Stats->Evaluate(base::t, base::U(), base::Work());
  

A full Problem.h example looks like this:

// -*- C++ -*-

#ifndef AMROC_PROBLEM_H
#define AMROC_PROBLEM_H

#define DIM  3
#define NEQUATIONS 15  // Euler equations for gases in 3D (5 fields), 
                       // +4 fields for passive scalars
                       // +1 field for temperature
                       // +1 field to output dcflag and 
                       // +1 field for subgrid kinetic energy 
                       // +2 field for subgrid stress components
                       // +1 field for sgs eps
#define NVARS     9    // Total number of active variables
#define NAUX      3

#define LESSVM_SGSEPS -2
#define LESSVM_SGS11 -8
#define LESSVM_SGS22 -7
#define LESSVM_SGS33 -6
#define LESSVM_SGS12 -5
#define LESSVM_SGS13 -4
#define LESSVM_SGS23 -3

#define CLES_IAUX_DEFINITION {LESSVM_SGS12, LESSVM_SGS23, LESSVM_SGSEPS}

#define OWN_AMRSOLVER

#include "WENOProblem.h"
#ifdef SYNCTIME
#include "WENOStdSyncTimeProblem.h"
#else
#include "WENOStdProblem.h"
#endif
#include "AMRInterpolation.h"
#include "WENOStatistics.h"

#define f_init_commongm FORTRAN_NAME(comblgm, COMBLGM)

extern "C" {
  void f_init_commongm(const INTEGER& meqn, DOUBLE* geom);
}

#ifdef OWN_AMRSOLVER

class SolverSpecific : 
#ifdef SYNCTIME
  public AMRSolverSyncTime<VectorType,FixupType,FlagType,DIM> {
  typedef AMRSolverSyncTime<VectorType,FixupType,FlagType,DIM> base;
#else
  public AMRSolver<VectorType,FixupType,FlagType,DIM> {
  typedef AMRSolver<VectorType,FixupType,FlagType,DIM> base;
#endif
  typedef AMRInterpolation<VectorType,DIM> interpolation_type;
  typedef F77FileOutput<VectorType,DIM> output_type;
  typedef WENOStatistics<VectorType,interpolation_type,output_type,DIM> stat_type;
public:
  SolverSpecific(IntegratorSpecific& integ, 
       base::initial_condition_type& init,
       base::boundary_conditions_type& bc) :
#ifdef SYNCTIME
    AMRSolverSyncTime<VectorType,FixupType,FlagType,DIM>(integ, init, bc) {
#else
    AMRSolver<VectorType,FixupType,FlagType,DIM>(integ, init, bc) {
#endif
    SetLevelTransfer(new F77LevelTransfer<VectorType,DIM>(f_prolong, f_restrict));
#ifdef f_flgout
    SetFileOutput(new F77FileOutput<VectorType,DIM>(f_flgout)); 
#else   
    SetFileOutput(new FileOutput<VectorType,DIM>()); 
#endif
    SetFixup(new FixupSpecific());
    SetFlagging(new FlaggingSpecific(*this)); 
    _Interpolation = new interpolation_type();
    _Stats = new stat_type(_Interpolation, (output_type*)_FileOutput);
  }  

  ~SolverSpecific() {
    delete _LevelTransfer;
    delete _Flagging;
    delete _Fixup;
    delete _FileOutput;
    delete _Interpolation;
    delete _Stats;
  }

  virtual void register_at(ControlDevice& Ctrl) { base::register_at(Ctrl); }
  virtual void register_at(ControlDevice& Ctrl, const std::string& prefix) {
    base::register_at(Ctrl, prefix);
    if (_Stats) _Stats->register_at(base::LocCtrl, ""); 
  }

  virtual void SetupData() {
    base::SetupData();
    f_init_commongm(DIM,geom);
    _Interpolation->SetupData(base::PGH(), base::NGhosts());
    _Stats->SetupData(base::PGH(), base::NGhosts(), base::shape, base::geom);
  }

  virtual void Advance(double& t, double& dt) {
        base::Advance(t, dt);
        _Stats->Evaluate(base::t, base::U(), base::Work());
    }

  virtual void Initialize_(const double& dt_start) { 
      base::Initialize_(dt_start);
      _Stats->Evaluate(base::t, base::U(), base::Work());
  }

protected:
  interpolation_type* _Interpolation; 
  stat_type* _Stats;
}; 
#endif

#endif

Sampling User Control

The user can control the behavior of the statistics package by specifing the following options within the AMRSolver section of solver.in:

Statistics {
      DefinitionFile stats.def
      OutputFile     stats
      Step           1000
   }

where the key DefinitionFile specifies the name of a text file containing all the definitions of statistcs to be extracted. OutputFile is the root name of the statistics files that will be produced. The timestep at which the statistics are written is attached to the root name. The key Step denotes the number of coarse steps to take between actual sampling.

Statistics Language Specification

The statistic language specified in DefinitionFile must conform to the following specification:

• each statistics ensamble is defined using the group keyword with matched curly braces enclosing all defintions.

       group {
              ....
       }
  

• you can specify an arbitrary number of comma separated expressions using the keys keyword.

       keys { exp1, exp2, exp3 }
  

• expresions are defined as arbitrary length mathematical expressions. The scanner/parser is able to recognize mixed Fortran and C style mathematical expressions. All variables accessible to the patch solver can be used (see bellow for a list). There are also a number of predefined constants that can be used.

       keys { u, v, w, u^2, v**2, u*w/2.0d0 }
  

• standard mathematic functions are also available and can be used within the expressions.
• an arbitrary number of one-dimensional samples of the predefined keys can the specified with the thread keyword. All threads are defined as a parametric curve with parametric coordinate a ranging from 0 to 1. The 2d or 3d coordinates of the points along the thread are sampled n times, at uniform intervals of a from 0 to 1. Two formats are accepted:

       thread { a, 0, 1 } @ 10
       thread { 0, a, 0 } @ <12>
  

  In the first case, 10 samplings are performed and the values are then reported in the output. In the second case, 12 samples are performed and the averaged such that a single value for each expression is written to the statistics file.
• analogous to threads, it is possible to define an arbitrary number of two-dimensional objects called surfaces, the corresponding keyword surface is then used. Like with threads, we must specify the bi-parametric representation of the surface in terms of a (the first parametric coordinate) and b (the second parametric coordinate). Like in the previous case, a and b vary from 0 to 1 in the number of steps specified, for example:

       surface {2*a, 1, b*zmax*39/40 } @ 20 <40>
       surface {0.845, 0.83+(1.17-0.83)*a, b*zmax*39/40 } @ <20> <40>
  

  In the first example, the second parametric coordinate is sampled over 40 points and the result is the averaged. So the output will contain only 20 points per expression. In the second case, both parametric coordinates are averaged, leading to a single value per expression.
• there is still a lightweight version of surface refered with the planes keyword. This allows a quick definition of planes aligned with the cartesian coordinates of the domain in a compact manner. Planes expect a comma separated definition of planes over which to perform statistics.

       planes  { x = xmin , x = 1, y = 0, z = [zmin , zmax/2] @ 20 }
       planes  { x @ nx, y @ 10 }  # equaly sampled planar statistics.
  

  In the first case, we define 23 planes, the first 3 are single plane definitions (arbitrary scalar expression can also be used here). The last 20 planes are defined as a range from zmin to zmax/2. All values samples on these planes are averaged so that the first plane definition is equivalent and will be executed literally as

       surface {xmin, ymin+(ymax-ymin)*a, zmin+(zmax-zmin)*b } @ <ny> <nz>
  

  In the second case, we use a similar notation but when the range is ommited, it is assumed that the range is from the minimum to the maximum of the domain size. This line also shows use of comments with the # keyword.

Here is a more general example of a statistics definition file

#
# Example of a statistics definition file. One should be able to get
# almost anything derived from the conservative or primitive vector
# of state with this approach.
#

group { # variables that will be measured.
        keys { u, v, w, rho * ((u^2 + v^2 + w^2)/2.0 + sgske), Y1, Y1^2,
               q(2)*q(3)/rho }
        # for 1d, 2d and 3d sampling (rack, surface and volume) we use
        # the dummy variables a and b (plane is a special form of surface).
        thread  { cos( 2 * pi * a), sin(2 * pi *a), 0 } @ 10
        surface { a, 1, b } @ 20 15
        surface { a, 2.,b } @ 20 <14>
        planes  { x = xmin , x = 1, y = 0, z = [zmin , zmax/2] @ 20 }
        planes  { x @ nx, y @ 10 }  # equaly sampled planar statistics.
        volume  # the whole volume.
        volume  { [0, 1], [0, 2], [1, 2] } # part of a volume
}

group {
        keys { T, p, E }
        surface {a + b, a - b, a } @ 20 20
}

the volume keyword is currently not implemented.

List of recognized symbols

-- CarlosPantano - 07 Dec 2005