Settings file

All interaction with Nemesis is done through settings files, which are in JSON-format. Note that whitespace is allowed and ignored around or between syntactic elements (values and punctuation, but not within a string value). Four specific characters are considered whitespace for this purpose: space, horizontal tab, line feed, and carriage return. The JSON syntax is derived from JavaScript object notation syntax, for which the following rules hold:

1. Data is in name/value pairs

2. Data is separated by commas

3. Curly braces hold objects

4. Square brackets hold arrays

Hint

The official JSON-format does not provide syntax for comments, but Nemesis allows comments and considers everything after ' // ' as a comment. These comments are stripped from the JSON string by the fileread_clean function, and that result is then passed on to the loadjson function (provided by the JSONLab module).

The main structure of the Nemesis settings file is given in the following listing:

{
  "nemesis" : {
  },
  "problem_data" : {
  },
  "function_data" : {
  },
  "mesh_data" : {
  },
  "field_data" : {
  },
  "physics_data" : {
  },
  "solver_data" : {
  },
  "output_data": {
  },
  "analysis" : {
  }
}

More details on the contents each of these objects can be found in their own sections.

After the settings file has been loaded by Nemesis (when the Nemesis is initialized, or by calling System.instance()), all settings can be accessed anywhere within Nemesis by using System.getSettings(). This function will return the entire structure with all settings. More specific settings can be accessed directly by passing an extra parameter to this function. For example: System.getSettings( 'problem_data.physical_parameters' ) will only return this part of the whole settings structure. This way it is also possible to get just a single return value instead of a structure.

Hint

If no settings file is passed, only the $nemesis_root/src/core/default_settings.json will be loaded. Passed settings files can overload these default settings. It is possible to initialize Nemesis by passing multiple settings files (or structures which correspond to the organization given here), where the order of them determines the final value for a specific setting (the last settings file / structure is the one used).

Block: nemesis

The nemesis block is used for Nemesis specific settings, see the following listing:

{
  "nemesis" : {
    "version"                           : 2.0,
    "verbose_level"                     : 1,
    "execution_mode"                    : "solve",
    "parallel_environment"              : "none",
    "pool_type"                         : "local",
    "pool_size"                         : 4,
    "pool_restart_forced"               : true,
    "tolerance"                         : 1e-12,
    "sparse_dense_limit"                : 0.05,
    "reuse_jacobians"                   : true,
    "show_warnings"                     : false
  }
}

The version indicates the Nemesis version, and it can be used to choose version specific methods during a run of Nemesis. The verbose_level is a positive (or zero) integer value, and it is commonly used to suppress output if set to 0. Larger values will (in general) provide more on screen output. The (optional) execution_mode defines what is done with the passed settings file. Valid options are "solve" (default), "init" (only load the settings file), or "post" (post-processing with already available data). The parallel_environment is currently not supported (only none is allowed), but it can be used in future versions to indicate to use the Matlab parallel pool (the pool_type, pool_size, and pool_restart_forced are necessary for this option) or the NMPI implementations.

Important

While MPI was supported in earlier versions (and by Keunsoo Park), during the implementation of unstructred meshes its support was lost. The overhead of the previous implementation was too large to be really efficient, while it did make everything else more complex to maintain. For a future re-implementation of parallel features one mainly needs to parallelize the FEHandler.m class. That class is used by a Field to assign the numbering of degrees of freedom (dofs) and quadrature points.

The tolerance is a Nemesis wide limit for nearly zero values (all values below this limit will be set to zero). The sparse_dense_limit (a value between 0 and 1) is used to determine when a dense array should be converted into a sparse array. Sparse arrays are only more efficient in their use if the number of nonzeros (nnz) is low. Here the limit is set to 5%. For efficiency, previously computed Jacobians can be reused or not (reuse_jacobians). Finally, show_warnings allows additional display of Matlab warning messages.

Block: problem_data

The problem_data block is used for problem specific settings, see the following listing:

{
  "problem_data" : {
   "problem_name"                      : "This is the name of the problem",
   "time" : {
     "method"                          : "BDF2",
     "start"                           : 0.0,
     "current"                         : 0.0,
     "end"                             : 2.0,
     "step_size"                       : 0.01,
     "step_number"                     : 0
   },

    "physical_properties" : {
      "sigma"                           : "sqrt(2) / 5"
    },
    "additional_properties" : {
      "use_kappa_LS_correction"         : true,
      "use_sharp_surface_tension"       : true,
      "use_sharp_density"               : true,
      "use_sharp_viscosity"             : true,
      "sharpening_method"               : "Heaviside",
      "interface_thickness"             : 0.125,
      "use_implicit_pressure"           : true
    },
    "verbose_level"                     : 0
  }
}

The problem_name is a reference that can be chosen freely by the user, it might be used for output of data. The time section determines if and how time is treated. Valid options for method are: "Steady" (default, time derivatives are set to zero), "Euler_Implicit", "BDF2" or "SpaceTime" (the mesh dimension will be increased with 1, so a 2D problem will create a 3D mesh). The time.start and time.end are the initial and final time levels, while the (automatically adjusted) time.current is the time value of the current time step. The time.step_size determines the number of steps required to reach the final time level. Finally, time.step_number is an automatically adjusted counter that indicates the (integer) number of steps taken to reach time.current. It should be equal to time.current / time.step_size.

Warning

If the time section is omitted the method will be set to Steady, so no advancement in time will be made.

Note

As BDF2 requires two previous time levels, for the first time step ImplicitEuler will be used. After this first time step BDF2 will be used again. Also the restart file will store all required previous solutions in time.

Attention

If a mesh is loaded from a Gmsh file, SpaceTime will add an extra dimension. The same mesh file can therefore be used for space-time and any of the time stepping methods.

The next two sections, physical_properties and additional_properties are quite important as they determine the nondimensional parameters used within the simulation. Possible options are

physical property	description	default value
length_scale	unit length scale	1
velocity_scale	unit velocity scale	1
time_scale	unit time scale	1
density	density value per phase	[1 1]
viscosity	viscosity value per phase	[1 1]
gravity	gravity vector	[0 0]
Cahn_number	Cahn number	0.06
Froude_number	Froude number	1
Peclet_number	Peclet number	100
Reynolds_number	Reynolds number	1
Weber_number	Weber number	1
sigma	surface tension	1
mobility	mobility (used by CH equation)	1
Weber_number	Weber number	1

Note

Some variables (like sigma probably) might also represent other physical properties, which will in general depend on the problem under consideration.

Block: function_data

The function_data block is used to define mathematical functions that can be used to initialize and set exact solution of Fields, see the following listing:

{
  "function_data" : {
    "droplet" : {
       "identifier" : 1,
       "class"      : "Interface",
       "parameters" : {
         "shape"  : "circle",
         "center" : [ 0.0, 0.0 ],
         "radius" : [ 1.0 ],
         "omega"  : "balanced"
       }
    },
    "exact" : {
      "identifier"  : 2,
      "class"       : "UserDefinedFunction",
      "parameters"  : {
        "u"   : "-sin(pi*x/2) .* cos(pi*y/2) .* cos(2*pi*t)",
        "v"   : " cos(pi*x/2) .* sin(pi*y/2) .* cos(2*pi*t)",
        "p"   : "zeros(x)",
        "rhs" : "exp( -( x.^2+y.^2 ) / $(problem_data.physical_properties.sigma)^2 )"
      }
    }
  }
}

Here the names droplet and exact are the function names, which can be used to as a reference to that specific function (by @function_name, see Block: field_data). The identifier is another way to reference the function, and should be a positive integer value. Next, class defines which Matlab m-file class (located in $(nemesis_root)/src/functions) must be used to process the parameters. The Interface class provides support for the variables c, omega, kappa and pressure, which are all used by the Navier-Stokes/Cahn-Hilliard method. For droplets circle/circular is used, and for ellipsoids ellipse/ellipsoid. Both require parameter values for the center and radius. Furthermore, flat interfaces require an additional parameters values for a point and normal vector.

Note

The number of elements of center, point and normal must correspond with the (spatial) mesh dimension. For ellipse shaped interfaces a radius value in each spatial direction is expected.

Furthermore, depending on the chosen Cahn-Hilliard method (see Block: physics_data), omega can be set. Valid options are classic (which leads to a non-uniform initial chemical potential, the traditional definition), balanced (which leads to an improved mass conservation, Kwakkel et al. (2020)), zero (uniform \(\omega = 0\)), one (uniform \(\omega = 1\)), or laplacianC (\(\omega = \Delta C\)).

More general functions can be constructed by the UserDefinedFunction class, where spatial coordinates (x, y, z) and time (t) can be used as variables. Furthermore, parameters defined elsewhere in the json file can be referenced by their position in the structure. For example, if a parameter sigma is defined under problem_data.physical_properties, it can be referenced by $(problem_data.physical_properties.sigma). Nemesis will automatically replace $(xyz) by the specific parameter xyz. The parameters of the UserDefinedFunction correspond to the variable provided by this function. For example, the given listing will provide u, v, p, and rhs. If the Field of variable u uses @exact, internally the function defined for u will be used: \(u = -sin(\pi x/2) \cdot cos(\pi y/2) \cdot cos(2\pi t)\) will be used.

Caution

It must be assured that the mathematical operators in the defined functions are valid Matlab operators. The variables x, y, z and t are arrays, and element-wise operators are therefore expected.

Block: mesh_data

The mesh_data block is used to define the geometry of the computational domain and the details of the finite elements used. An example is given in the following listing:

{
   "mesh_data" : {
    "geometries" : [{
      "identifier"                      : 1,
      "type"                            : "octree",
      "maximum_level"                   : 1,
      "refinement_ratio"                : [ 2, 2 ],
      "refinement_conditions"           : [ { "field"      : "c",
                                              "solution"   : "exact",
                                              "range"      : [ 0.1, 0.9 ],
                                              "action"     : [-1,1,-1]
                                            },
                                            { "coordinate" : "x",
                                              "range"      : [ 0.9, 1.1 ],
                                              "action"     : [-1,1,-1]
                                             }
                                          ],
      "mesh_refinement"                 : {
        "verbose_level"           : 0,
        "prevent_hanging_nodes"   : false
      },
      "input_type"                      : "matlab",
      "input_filename"                  : "unstructured_square_2x2_L1.msh",
      "mesh_name"                       : "channel",
      "number_of_elements"              : [  8, 8 ],
      "x"                               : [ -2, 2 ],
      "y"                               : [ -2, 2 ],
      "periodic"                        : [ false, false ],
      "patch_names"                     : [ "x0y0", "x1y0", "x0y1", "x1y1",
                                            "left", "right", "bottom", "top",
                                            "fluid" ],
      "verbose_level"                   : 0
    }],
    "finite_elements" : [{
      "identifier"                      : 1,
      "polynomial_type"                 : [ "hermite", "hermite" ],
      "continuity"                      : [ 1, 1 ],
      "number_of_modes"                 : [ 5, 5 ],
      "quadrature_order"                : [ 6, 6 ],
      "refinement_ratio"                : [ 2, 2 ]
    }],
    "time_element" : {
      "polynomial_type"                 : "Lagrange",
      "continuity"                      : 0,
      "number_of_modes"                 : 4,
      "quadrature_order"                : 6,
      "refinement_ratio"                : 1
    }
  }
}

Note

The square brackets make it possible to define multiple geometries, where each geometry must be enclosed within curly brackets. The same holds for the finite_elements, where each (spatial) finite element must be enclosed within curly brackets.

The identifier is a unique reference number of the geometry (it might be used to identify a specific geometry). The type must be "static" for static meshes (no adaptive mesh refinement), or "octree" if octree-based adaptive mesh refinement is possible. If the latter is used, one must also define maximum_level, which is an postive (or zero) integer value indicating how many refinement levels are allowed.

Note

If the maximum_level is set to zero, no adaptive mesh refinement is used. However, the mesh object will be constructed from the OctreeMesh class anyway. In principle StaticMesh and OctreeMesh (with maximum_level = 0) must given identical results.

The refinement_ratio sets the number of child elements in each spatial direction. Default is a refinement factor of 2 in each spatial direction (which leads to a quadtree in 2D and an octree in 3D meshes).

Warning

The refinement_ratio here must correspond with the refinement_ratio of the finite element!

The refinement_conditions provide the rules where adaptive mesh refinement will occur. The regions of the mesh that meet the set conditions will be refined up to maximum_level. It is possible to define multiple rules (by using square and curly brackets), and each rule can be based on a field or on mesh coordinates. In case of a field, valid options to pass are: field with the variable name of the field, solution to specify a specific solution of the field (if not specified last will be used), range to set the lower and upper limits for the refinement condition, and action which specifies what should happen (-1 for coarsening, +1 for refinement) for elements where element_value < lower_limit (if refined), lower_limit < element_value < upper_limit, (if the element refinement level < maximum_level) or element_value > upper_limit (if refined). The default is [-1,1,-1], which means that only elements within the range will be refined and will be coarsened again once they are outside this range.

Tip

By default the field value will be used, but the variable name of the field can also be used to use gradients as conditions. The specific gradient can be set by including it in the variable name as follows u.x (for the x-gradient of the u variable), or p.yy (for the second derivative in y-direction of the p variable). To be able to use a variable is must be defined as a field, and to use second order gradients the underlying finite elements must support that (for second order derivatives Hermite basis function must be used).

Warning

Other values than [-1,1,-1] have not been tested properly, so be careful when using other settings.

Attention

In a future version all options related with the refinement could/should be positioned within the mesh_refinement block.

The prevent_hanging_nodes option is still an alpha option, and used to apply a special refinement that splits an element into 3 child elements.

From which source the mesh is created by the input_type variable, and valid options are matlab or gmsh. For the latter an input_filename must be provided, which should be a valid Gmsh file with the msh extension. In this case all Matlab mesh options are not required as the mesh file will provide all necessary information. If the matlab mesh is used, the following options are possible:

Matlab mesh options	description	required?
mesh_name	channel or circular	yes
equidistant	GLL element layout if false	yes
x, y, z	mesh limits in each spatial directions	yes (up to mesh dimension)
number_of_elements	mesh dimensions in each spatial direction	yes (up to mesh dimension)
periodic	mesh periodicity in each spatial directions	yes (up to mesh dimension)
patch_names	overload the internally created patch_names	no

Block: field_data

The field_data block is used to define the variables (Field objects in Matlab) and their soltuions (Solution objects in Matlab) which are represented on a mesh and are used by the equations (see Block: physics_data) to set previous values and gradients. By solving a set of equations, the solutions in the corresponding fields will be updated. The following listing provides an example of the definition of two variables: u and residual.

{
   "field_data" : {
     "verbose_level" : 0,
     "initialize"    : true,
     "output"        : [ "u", "residual" ],
     "fields"        : [
       { "variable"                        : "u",
         "mesh_identifier"                 : 1,
         "fe_handler_data"                 : {
           "identifiers"             : 1,
           "shared_dofs"             : "all",
           "shared_quadrature_nodes" : false,
           "shared_finite_elements"  : true
         },
         "solution_data"                   : [ { "type" : "initial", "value" : 0 },
                                               { "type" : "exact"  , "value" : "@exact" },
                                               { "type" : "error"  , "value" : 0, "init_with_bc": false } ],
         "boundary_data"                   : {
           "fluid"     : { "type":"Internal" , "value": ["xx","yy"], "weight" : [100,100] },
           "left"      : { "type":"Dirichlet", "value":"@exact", "weight" : 100 },
           "right"     : { "type":"Dirichlet", "value":"@exact", "weight" : 100 },
           "top"       : { "type":"Dirichlet", "value":"@exact", "weight" : 100 },
           "bottom"    : { "type":"Dirichlet", "value":"@exact", "weight" : 100 }
         },
         "monitor_data" : { "filename"     : "u_data.dat",
                            "write_header" : true,
                            "variables"    : [ "timestep", "time", "error.norms.L2", "number_of_dofs" ] }
       },
       { "variable"                        : "residual",
         "mesh_identifier"                 : 1,
         "fe_handler_data"                 : {
           "identifiers"             : 1,
           "shared_dofs"             : "none",
           "shared_quadrature_nodes" : false,
           "shared_finite_elements"  : true
         },
         "solution_data"                   : [ { "type" : "initial", "value" : 0 } ]
       }
     ]
   }
}

Here verbose_level allows to select a specific (debug) verbosity during the use of the fields. Next, initialize defines if the fields need to be initialized (default is true). The variable names defined by output will be exported to a Paraview VTK file during runtime. The most important part of the field_data block is the fields block, which is an array (square brackets) of variables to be defined. Each independent field is defined between curly brackets. Detailed specifications of fields, solution_data, boundary_data, and monitor_data are given below.

fields specification

A Field object represent a variable, and internally it is defined on a mesh (depends on a mesh object) and one (or more) finite elements (finite element object(s)). A FiniteElement is the mapping between physical and natural domain (standard element domains ranging between -1 and +1). By default, a single FiniteElement is used within each mesh_element, but any support for multiple finite elements (necessary for p-refinement) can be added (through the FEHandler object within the Field). Practically, a variable can have multiple solutions of interest. Therefore, within the Field object multiple Solution objects are created which can represent the last, exact, error, nonlinear, coupling, last-n. The last solution is the most recent solution, and at initialization initial will be converted to last automatically. If there is an exact solution, also an error solution will be defined, which is the difference between last and exact. The nonlinear/coupling solutions are temporary solutions, used to check if nonlinear/coupling convergence has occurred. Note that these temporary solution will be deleted again after the solution has convergence. Finally, last-n represents solutions at time level last-n, which can be used for time stepping purposes.

Note

It might be more logical/efficient to include the location of the quadrature points to the mesh object as well, as these locations are only present in the physical domain. All Jacobian operations could/should then be moved to the mesh object as well. The Field object can then be defined by a polynomial type (Lagrange, Hermite) and polynomial order. A simple change of the reference mesh (quadrature order / number of quadrature points) could then be used to represent/plot a solution on a more detailed mesh. The finite elements still form the ‘bridge’ between the mesh (physical domain) and the basis function weights (alpha values, stored in the Field/Solution).

Field options	description	required?	remarks
variable	Name of the variable	yes
mesh_identifier	ID of the underlying Mesh object	yes
periodic	Follow the periodicity of the Mesh object	no (default is true)	if the Mesh object is not periodic, this won’t do anything
fe_handler_data	Detailed information on the FiniteElements	yes	see the fe_handler_data specification
solution_data	Block with Solutions (see below)	yes (at least initial)	see the solution_data specification
boundary_data	Boundary conditions of this Field	no	see the boundary_data specification
monitor_data	Values to export for this Field	no	see the monitor_data specification

fe_handler_data specification

Each Field has an FEHandler object, which controls the assignment of finite elements, performs numbering of degrees of freedom and quadrature nodes and constructs all operators (like H, Dx, Dxx, etc., and also the internal boundary operators). The options for the FEHandler object of a field are defined in fe_handler_data, and its specification is as follows:

fe_handler_data options	description	required?	remarks
identifiers	ID(s) of the FiniteElements	yes	single ID or array with IDs
shared_dofs	Which DOFs should be shared?	no (default is "all")	valid options: "all", "values", "none"
shared_quadrature_nodes	Should quadrature nodes be shared?	no (default is true)	see the second note below
shared_finite_elements	Should finite elements be shared?	no (default is true)	reuse a single, or make a cell array of FiniteElements

Note

By default all DOFs (numbering) will be shared between elements, which means for C0 basis functions that the values will be shared and for C1 basis functions that both the values and first derivatives will be shared (and the mixed derivatives). The values option will limit the sharing of DOFs for C1 basis functions to only the DOF controlling the value. Any derivative DOFs that would be shared with all are now decoupled between elements. Finally, none will ensure no DOFs are shared between elements. For both reduced DOF sharing methods the coupling between the elements must be ensured through internal boundary conditions, see boundary_data specification. Without proper internal boundary conditions the solver will probably complain with rank deficient warnings or just not work at all. If the warnings occur even with proper internal boundary conditions, the quadrature order might be too low.

Note

Consider C0 Lagrange basis functions defined on two 1D (edge) elements A and B. At their mutual vertex (point C), only the value will be continuous. This C0 continuity is achieved in a strong sense, as the ID of the DOF is shared between the two elements. The first derivative at point C is generally not continuous, as it can be defined as a combination of the basis functions within element A or as a combination of the basis functions within element B. However, one could enforce continuity by combining the basis functions of both elements to define a single derivative value at point C. The latter approach is enabled in the json file by setting shared_quadrature_nodes to true.

solution_data specification

A Field can contain multiple solutions (Solution objects), and the specification of each solution is as follows:

Solution options	description	required?	remarks
type	Type of the solution (string)	yes	initial will be converted automatically to last
value	Value to assign at initialization	yes	It is possible to make a reference to a UserDefinedFunction
init_with_bc	Apply boundary conditions at initialization?	no (default is true)

boundary_data specification

A Field object can have its own boundary conditions, which only depend on the field itself and are the same for all (internal) solutions (although init_with_bc will allow solutions without a direct connection to the boundary conditions). The specification of the boundary_data is as follows:

boundary_data options	description	required?	remarks
patch_name	Patch name to for the condition (variable)	yes	The patch should be available for the mesh object
type	Boundary condition type	yes	Dirichlet, Neumann (see FEBoundary.m)
value	Value for the boundary condition	yes	It possible to make a reference to a UserDefinedFunction
form	Weak or strong boundary condition	no (default is weak)	Strong boundary conditions are not supported currently
weight	Weight in the LS functional (for weak bc only)	yes	Zero weight will disable the boundary condition

In order to make a reference a UserDefinedFunction, that function must provide a parameter with the variable name. For the given example here, the exact function must have a parameter with the name u, see the Block: function_data. A special boundary condition type is the internal boundary condition, which only has been used so far for 2D meshes (any other use might result in errors). The patch_name for the internal boundary condition is the Patch containing all 2D elements, in the example given here fluid is used. The type is set to Internal, which will ensure the use of the setInternal function in the FEBoundary object. The value is a block of strings with the operators to make continuous between elements, valid strings here are: value, x, y, xx, xy, yy (given the chosen polynomial is able to provide these). Finally, weight contains a block of numeric values with the least-squares weights for each provided operator (the number of weights must exactly match the number of operators!). In the given example, the xx-derivative, and yy-derivative between elements is enforced in a weak sense (each has a weight of 100).

Warning

For the internal boundary condition there is also a periodic option, which is set to false by default. Although it might work, enabling periodicity has not been verified.

monitor_data specification

Output of a singular Field is controlled through the monitor_data options. These are given by:

monitor_data options	description	required?	remarks
enabled	Boolean to enable / disable the monitor	no	By default the monitor is enabled (true)
filename	Filename to write to	no	It will be written to the general_data directory (default: $variable_data)
write_header	Boolean to enable / disable the header	no	By default the header will be written (true)
header	Add a user defined header	no	If not provided, the header will contain the variable names
format	Write variables with a specific format	no	Matlab formatSpec is expected (for example: “%5d” or “%8.3f”)
variables	Block with variables to write to the file	yes	These variables should be available in the object

Note

If a variable is not available in the Field object, an error will terminate the execution of Nemesis. Maybe add some basic try/catch support to the Monitor class (it is a superclass of the Field class).

Block: physics_data

The physics_data block is used to define the physics (equations) to be solved by Nemesis. These physics can be coupled, but this is not a requirement. The following listing provides an example of the definition of three systems (equations). The top level groups are systems, coupling, verbose_level. The first defines the actual physics/equations to solve for, which will be explained in detail later. The coupling group sets the maximum number of coupling iterations, and an option could be added (not supported at the moment) that sets which Fields are used to determine coupling convergence. The verbose_level is a flexible option to determine the verbosity of the output, but no predefined output options or guidelines have been implemented for this at the moment (general rule is that a lower value results in less output).

{
  "physics_data" : {
    "systems" : [
    { "class_name"                        : "NavierStokes",
      "identifier"                        : 1,
      "variables"                         : ["u" "v" "p"],
      "coupled_variables"                 : ["c", "curvature" ],
      "number_of_equations"               : 3,
      "equation_weights"                  : [ 0.01, 0.01, 1 ],
      "nonlinear" : {
        "linearization" : "Newton",
        "max_iteration" : 5
      },
      "boundary_data"                   : {
        "top"       : { "type"            :"gnbc",
                        "value"           : 0,
                        "interface_field" : "c",
                        "velocity"        : ["u","v"],
                        "beta_f"          : "$(friction_coefficient)",
                        "Uw"              : "$(wall_velocity)",
                        "theta_0"         : "$(equilibrium_contact_angle)",
                        "Cn"              : "$(Cahn_number)",
                        "Ca"              : "$(Capillary_number)",
                        "alpha"           : "$(alpha)",
                        "gprime"          : "@switching_function",
                        "weight"          : 100
                      },
        "bottom"    : { "type"            :"gnbc",
                        "value"           : 0,
                        "interface_field" : "c",
                        "velocity"        : ["u","v"],
                        "beta_f"          : "$(friction_coefficient)",
                        "Uw"              : "-$(wall_velocity)",
                        "theta_0"         : "$(equilibrium_contact_angle2)",
                        "Cn"              : "$(Cahn_number)",
                        "Ca"              : "$(Capillary_number)",
                        "alpha"           : "$(alpha)",
                        "gprime"          : "@switching_function",
                        "weight"          : 100
                      }
      },
      "residual_fields"                   : [ "residual_NS" ],
      "properties"                        : {
        "pressure_gradient"     : [ 0, 0 ],
        "epsilon_div_pressure"  : 0.00
      },
      "mesh_identifier"                 : 1,
      "verbose_level"                   : 0
    },
    { "class_name"                        : "CahnHilliard",
      "short_name"                        : "CH",
      "identifier"                        : 2,
      "potential_type"                    : "balanced",
      "variables"                         : [ "c", "omega", "kappa" ],
      "coupled_variables"                 : [ "u", "v" ],
      "residual_fields"                   : [ "residual_CH" ],
      "number_of_equations"               : 2,
      "equation_weights"                  : [ 1, 1 ],
      "mobility_type"                     : 1,
      "explicit_divergence"               : false,
      "nonlinear" : {
        "linearization" : "Newton",
        "max_iteration" : 5
      },
      "mesh_identifier"                 : 1,
      "verbose_level"                   : 0
    },
    { "class_name"                        : "EnergyEquation",
      "short_name"                        : "EE",
      "identifier"                        : 3,
      "variables"                         : ["T"],
      "coupled_variables"                 : ["u","v","c"],
      "number_of_equations"               : 1,
      "equation_weights"                  : [ 1 ],
      "residual_fields"                   : [ "residual_E" ],
      "mesh_identifier"                   : 1,
      "verbose_level"                     : 0
    }],
    "coupling" : {
      "max_iteration" : 5
    },
    "verbose_level" : 0
  },
}

The systems group is a block array in which each physics/equation is added within curly brackets. Each of these has a similar structure, so only the first one (NavierStokes) will be explained here in detail. The

keyword	description
class_name	Name of the m-file that controls/defines this physics (the spatial dimension of the mesh will be used to select equation dimension)
variables	Block array with variable names to solve for (must be defined in field_data)
coupled_variables	Block array with variable names to use for coupling (must be defined in field_data)
number_of_equations	Number of equations defining the physics (it allows to select different combinations of equations)
equation_weights	Block array with weights for each equation (size must match number_of_equations)
nonlinear	Defines the settings for nonlinear equations (linearization method and max_iteration)
boundary_data	Boundary data for the physics, see the note below.
residual_fields	Fields to use to store the residual of this physics (must be defined in field_data)
properties	Additional properties required for this physics only (optional, depends on the physics)
mesh_identifier	Identifier of the Mesh object to solve this physics on (see the note below)
verbose_level	Control the verbosity of the output of this physics

Important

In general boundary_data (boundary conditions) will be set on the Fields directly (see Block: field_data for more details how). However, if a boundary condition depends on multiple Fields (like the Generalized Navier Boundary Condition, GNBC) this can’t be used. In that case the boundary condition (for the Patch it works on) must removed from the Fields, and be replaced by a boundary condition for the Physics. These boundary conditions are added under boundary_data. The general structure of a boundary conditions remains the same as used for a Field, but more detailed parameters can be passed to the boundary condition.

Note

Internally Nemesis will call a method within the Physics class with the name bc_TYPE, where TYPE is the type value given in the structure. In the given example from above: the type is set to gnbc. Since this boundary condition is defined within the NavierStokes group, Nemesis will call the NavierStokes.bc_gnbc( patch_name, bc_data ) method. This method is the actual implementation of the Generalized Navier Boundary Condition (GNBC). The expected patch_name is the corresponding patch name (boundary name) the condition works on (here top and bottom), and bc_data is the full structure as defined as the value of the patch name. Here bc_data will be the structure defined by boundary_data.top and boundary_data.bottom). This design makes it quite easy to access anything from the JSON in the actual boundary condition.

Note

For most operations the identifier variables are not really used (only a single identifier is supported). However, theoretically it is possible to solve different physics on different meshes (which might lead to interpolation). Or one could use different finite elements for different meshes, or different finite elements within a single mesh. Although most of these specialized options are not supported at the moment, identifiers might come in handy if someone is will to invest time in the implementations.

Block: solver_data

The solver_data block is used to define the specifics of the solver that will be used. This block will be passed to the constructor of the LeastSquaresSolver class, or when the solve method of that class is called.

{
  "solver_data" : {
    "solver_method"           : "direct",
    "reorder_method"          : "none",
    "direct_settings" : {
      "least_squares_method"  : "QR"
    },
    "iterative_settings" : {
      "method"                : "pcg",
      "use_mex"               : false,
      "convergence_criterion" : 1e-12,
      "max_iteration"         : 1000,
      "preconditioner"        : {
        "method"    : "ichol",
        "type"      : "nofill",
        "michol"    : "off",
        "droptol"   : 1e-2,
        "diagcomp"  : 1.0,
        "shape"     : "lower"
      },
      "verbose_level"         : 0
    },
    "remove_zero_columns"     : true,
    "verbose_level"           : 0
  }
}

The solver_method can be set to direct or iterative. If direct is used, the direct_settings.least_squares_method. can be set to NE or QR. Here NE indicates that the normal equations will be created, and QR indicates that so-called direct minimization will be used (solved by a QR decomposition of the matrix). For more details on the two approaches, see [8]. The reorder_method might speed up the solving of the system of equations. For more details on the valid options (none, colamd, symamd, symrcm) see Matlab reorder options. If an iterative solver is selected, the iterative_settings must contain all settings of the iterative solver. Typical properties to set are method (default is pcg, the Preconditioned Conjugate Gradients method), convergence_criterion, the maximum number of iterations (max_iterations), and the preconditioner to use. For the latter, see the documentation of Matlab on iterative methods for linear systems. The possible options for the preconditioner match the ones available by default in Matlab. In general, verbose_level allows the user to set the verbosity of the solver.

Note

The remove_zero_columns option will explicitly remove any columns in the matrix that contain all zeros (only use when the solver does complain it needs this setting).

Block: output_data

The output_data block is used to set the specifics for output of data during the simulation, see the following listing:

{
  "output_data": {
    "verbose_level"             : 1,
    "enable_output"             : true,
    "root_directory"            : "./two_phase_heat_transfer/",
    "general_data" : {
      "directory"               : "./post",
      "write_session"           : false,
      "session_filename"        : "session.out",
      "max_number_of_backups"   : 2
    },
    "solution_data" : {
      "format" : "paraview",
      "directory" : "./vtk_data",
      "filename" : "solution_$(problem_data.problem_name)_mesh$(mesh_id)",
      "write_frequency" : 1
    },
    "restart_data" : {
      "directory"           : "./restart",
      "filename"            : "restart",
      "write_frequency"     : 1
    }
  }
}

The enable_output is the main switch to control any output. In order to generate any output, it needs to be set to true. The root_directory is the top-level output location, and it can be used as a relative location with respect to the json file. Next there are 3 groups of data: general_data, solution_data, and restart_data. Of these, general_data will contain all output like the Matlab session file and any output from monitors. In general, an existing session file will not be overwritten by repeated use of the same json file. Each new run will backup any existing session files with a number. To limit the number of backups, max_number_of_backups can be used. The second, solution_data, will contain output of the solutions, which should be in the Paraview as shown. Other data format options can be added by the user of course. The write_frequency value will determine if output is requested every time step (value==1) or if output should be less frequent. Finally, the restart_data is used to write data files which can be used to restart a simulation. By default the restart data is included in the json files written to the restart_data.directory. All required data of fields and solutions is added as binary data within the json. This allows the restart data to be loaded from a single file.

Block: analysis

The analysis_data block can be used for analysis purposes. However, up to now it has not been used as such. Therefore, it can be left out for regular users.