PODReducedBasisTrainer

Computes the reduced subspace plus the reduced operators for POD-RB surrogate.

Overview

In this trainer, an intrusive Proper Orthogonal Decomposition (POD) based Reduced Basis (RB) method (Pinnau (2008)) is implemented which, unlike non-intrusive surrogates such as Polynomial Regression or Polynomial Chaos Expansion, is capable of considering the physics of the full-order problem at surrogate level. Therefore, it is often referred to as a physics-based but still data-driven approach. The intrusiveness, however, decreases the range of problems which this method can be used for. In the current version, this surrogate model can deal with Parameterized scalar-valued linear steady-state PDEs with affine parameter dependence only. This class is responsible for two steps in the generation of the surrogate model:

  1. Generation of reduced subspaces

  2. Generation of reduced operators

It must be mentioned that in POD-RB literature, the training phase is often referred to as offline phase and in the upcoming sections the two expressions are used interchangeably.

A parameterized linear steady-state PDE

Before the details of the above-mentioned steps are discussed, a short overview is given about the problems considered. A scalar-valued linear steady-state PDE can be expressed in operator notation as:

Au=b,\mathcal{A}u = \mathcal{b},(1)

where uu is the solution, A\mathcal{A} is a linear operator and b\mathcal{b} is a source term. The linear operator and the source terms may depend on uncertain parameters which are denoted by μi, i=0,...,Nμ\mu_i,~i=0,...,N_\mu and organized into a parameter vector μ=[μ1,...,μNμ]T\boldsymbol{\mu} = [\mu_1,...,\mu_{N_\mu}]^T. Therefore, Eq. (1) can be expressed as:

A(μ)u=b(μ).\mathcal{A}(\boldsymbol{\mu})u = \mathcal{b}(\boldsymbol{\mu}).(2)

This also means that the solution itself is the function of these parameters u=u(μ)u=u(\boldsymbol{\mu}). To make an efficient surrogate, operator A(μ)\mathcal{A}(\boldsymbol{\mu}) and source b(μ)\mathcal{b}(\boldsymbol{\mu}) should have an affine parameter dependence:

A(μ)=i=1NAfiA(μ)Ai,and b(μ)=i=1Nbfib(μ)bi,\mathcal{A}(\boldsymbol{\mu}) = \sum \limits_{i=1}^{N_A} f^A_i(\boldsymbol{\mu})\mathcal{A}_i, \text{\quad and \quad} \mathcal{b}(\boldsymbol{\mu}) = \sum \limits_{i=1}^{N_b} f^b_i(\boldsymbol{\mu})\mathcal{b}_i,

or in other words, the operators have to be decomposable as the sums of products of parameter-dependent scalar functions and parameter-independent constituent operators. By plugging the decompositions back to Eq. (2), the problem takes the following form:

(i=1NAfiA(μ)Ai)u=i=1Nbfib(μ)bi.\left(\sum \limits_{i=1}^{N_A} f^A_i(\boldsymbol{\mu})\mathcal{A}_i\right)u = \sum \limits_{i=1}^{N_b} f^b_i(\boldsymbol{\mu})\mathcal{b}_i.(3)

Therefore, before starting the construction of a POD-RB surrogate model, the user must identify the these decompositions first. At input level, this can be done by utilizing the Tagging System. For each constituent operator a separate vector tag has to be created and the tags need to be supplied to the trainer object through the "tag_names" input parameter. Furthermore, an indicator shall be added to each tag through the "tag_types" input to show if the tag corresponds to a source term (bi\mathcal{b}_i) or an operator (Ai\mathcal{A}_i). As a last step, this system is discretized in space using the finite element method to obtain:

(i=1NAfiA(μ)Ai)u=i=1Nbfib(μ)bi,\left(\sum \limits_{i=1}^{N_A} f^A_i(\boldsymbol{\mu})\textbf{A}_i\right)\textbf{u} = \sum \limits_{i=1}^{N_b} f^b_i(\boldsymbol{\mu})\textbf{b}_i,(4)

where Ai\textbf{A}_i and bi\textbf{b}_i are finite element matrices and vectors, while u\textbf{u} denotes the vector containing the values of the degrees of freedom. This model is referred to as Full-Order Model (FOM) in subsequent sections. Furthermore, let u(μ)\textbf{u}(\boldsymbol{\mu}^* ) denote a solution vector which is obtained by solving Eq. (4) with μ=μ \boldsymbol{\mu}=\boldsymbol{\mu}^* ~.

Generation of reduced subspaces

Even though it is not explicitly stated in the previous sub-section, u\textbf{u} may contain solutions for multiple variables, hence it can be expressed as u=[u1;...; uNv]\textbf{u}=[\textbf{u}_1;...;~\textbf{u}_{N_v}], where NvN_v is the total number of variables. It is assumed that each variable has Ni, i=1,...,NvN_i,~i=1,...,N_v spatial degrees of freedom, thus the size of the full-order system is i=1NvNi\sum\limits_{i=1}^{N_v}N_i.

As a first step in this process, Eq. (4) is solved using NsN_s different parameter samples and the solution vectors for each variable defined on the "var_names" input parameter are saved into snapshot matrices

Si=[ui(μ1),...,ui(μNs)].\textbf{S}_i = [\textbf{u}_i(\boldsymbol{\mu}_1),...,\textbf{u}_i(\boldsymbol{\mu}_{N_s})].(5)

In this implementation, the solutions are obtained from a PODFullSolveMultiApp using a PODSamplerSolutionTransfer and the snapshot matrices are stored within the trainer object. The next step in this process is to use these snapshots to create reduced sub-spaces for each variable. This can be done by performing POD on the snapshot matrices, which consists of the following four steps for each variable:

  1. Creation of correlation matrices: The correlation matrices (Ci\textbf{C}_i) can be computed using the snapshot matrices as Ci=SiTWiSi\textbf{C}_i=\textbf{S}_i^T \textbf{W}_i \textbf{S}_i, where Wi\textbf{W}_i is a weighting matrix. At this moment only Wi=I\textbf{W}_i = \textbf{I} is supported, where I\textbf{I} is the identity matrix.

  2. Eigenvalue decomposition of the correlation matrices: The eigenvalue decompositions of the correlation matrices is obtained as: Ci=ViΛiViT\textbf{C}_i = \textbf{V}_i\boldsymbol{\Lambda}_i\textbf{V}^T_i, where matrix Vi\textbf{V}_i and matrix Λi\boldsymbol{\Lambda}_i contain the eigenvectors and eigenvalues of Ci\textbf{C}_i, respectively.

  3. Determining the dimension of the reduced subspace: Based on the magnitude of the eigenvalues (λi,k, k=1,...,Ns\lambda_{i,k},~k=1,...,N_s) in Λi\boldsymbol{\Lambda}_i, one can compute how many basis functions are needed to reconstruct the snapshots with a given accuracy. The rank of the subspace (rir_i) can be determined as:

    ri=arg min1riNs(k=1riλi,kk=1Nsλi,k>1τi),r_i = \argmin\limits_{1\leq r_i\leq N_s} \left(\frac{\sum_{k=1}^{r_i}\lambda_{i,k}}{\sum_{k=1}^{N_s}\lambda_{i,k}} > 1 - \tau_i\right),(12)

    where τ\tau is a given parameter describing the allowed error in the reconstruction of the snapshots. This can be supplied to the trainer using the "error_res" input parameter. Of course, this rank can be determined manually as well. For more information about this option, see PODReducedBasisSurrogate.

  4. The reconstruction of the basis vectors for each variable: For this, the eigenvalues and eigenvectors of the correlation matrices are used together with the snapshots as:

    Φi,k=1λi,kj=1NsVi,k,jSi,j,\Phi_{i,k}=\frac{1}{\sqrt{\lambda_{i,k}}}\sum\limits_{j=1}^{N_s}\textbf{V}_{i,k,j}S_{i,j},(13)

    where Φi,k\Phi_{i,k} is the kk-th (k=1,...,rik=1,...,r_i) basis function of the reduced subspace for variable ui\textbf{u}_i. Moreover, Vi,k,j\textbf{V}_{i,k,j} denotes the jj-th element of the kk-th eigenvector of correlation matrix Ci\textbf{C}_i. It is important to remember that Φi,k\Phi_{i,k} has a global support in space, and shall not be mistaken for the local basis functions (ϕ\phi) of the finite element approximation. The global basis vectors can be also referred to as POD modes and the two expressions are used interchangeably from here on.

It must be noted that these basis functions are also stored in the trainer object. Finally, the solutions of different variables in the full-order model can be approximated as the parameter-dependent linear combination of these basis functions:

ui(μ)k=1riΦi,kci,k(μ),or ui(μ)Φici(μ)\textbf{u}_i(\boldsymbol{\mu})\approx\sum\limits_{k=1}^{r_i}\Phi_{i,k}c_{i,k}(\boldsymbol{\mu}), \text{\quad or \quad} \textbf{u}_i(\boldsymbol{\mu})\approx\boldsymbol{\Phi}_{i}\textbf{c}_{i}(\boldsymbol{\mu})(6)

where Φi=[Φi,1,...,Φi,ri]\boldsymbol{\Phi}_{i}=[\Phi_{i,1},...,\Phi_{i,r_i}] is a matrix with the POD modes as columns and ci(μ)=[ci,1(μ),...,ci,ri(μ)]T\textbf{c}_{i}(\boldsymbol{\mu}) = [c_{i,1}(\boldsymbol{\mu}),...,c_{i,r_i}(\boldsymbol{\mu})]^T are the expansion coefficients. In essence, these coefficients describe the coordinates of the approximate solution in the reduced subspace. To approximate the full solution vector u(μ)\textbf{u}(\boldsymbol{\mu}) using its components (ui(μ)\textbf{u}_i(\boldsymbol{\mu})-s), a segregated approach is used as follows:

u(μ)=[u1uNv][Φ1c1(μ)ΦNscNs(μ)]=[Φ1000Φ2000ΦNs][c1(μ)c2(μ)cNs(μ)]=Φc(μ).\textbf{u}(\boldsymbol{\mu}) = \left[ \begin{array}{c} \textbf{u}_1\\ \vdots \\ \textbf{u}_{N_v} \end{array}\right] \approx \left[ \begin{array}{c} \boldsymbol{\Phi}_{1}\textbf{c}_{1}(\boldsymbol{\mu})\\ \vdots \\ \boldsymbol{\Phi}_{N_s}\textbf{c}_{N_s}(\boldsymbol{\mu}) \end{array}\right] = \left[ \begin{array}{cccc} \boldsymbol{\Phi}_{1} & 0 & \cdots & 0\\ 0 & \boldsymbol{\Phi}_{2} & \cdots & 0\\ \vdots & \vdots & \ddots & \vdots\\ 0 & 0 & \cdots & \boldsymbol{\Phi}_{N_s} \end{array}\right] \left[ \begin{array}{c} \textbf{c}_{1}(\boldsymbol{\mu})\\ \textbf{c}_{2}(\boldsymbol{\mu})\\ \vdots\\ \textbf{c}_{N_s}(\boldsymbol{\mu}) \end{array}\right] = \boldsymbol{\Phi}\textbf{c}(\boldsymbol{\mu}).(7)

It is important to mention that in this approximation the unknowns are the elements of c(μ)\textbf{c}(\boldsymbol{\mu}) vector. In most of the cases, the size of this vector (i=1Nvri\sum\limits_{i=1}^{N_v}r_i) is considerably smaller than the size of the original solution vector (i=1NvNi\sum\limits_{i=1}^{N_v}N_i).

Generation of reduced operators

To generate reduced opertors, Eq. (7) is plugged into Eq. (4) first:

(i=1NAfiA(μ)Ai)Φc(μ)=i=1Nbfib(μ)bi.\left(\sum \limits_{i=1}^{N_A} f^A_i(\boldsymbol{\mu})\textbf{A}_i\right)\boldsymbol{\Phi}\textbf{c}(\boldsymbol{\mu}) = \sum \limits_{i=1}^{N_b} f^b_i(\boldsymbol{\mu})\textbf{b}_i.(8)

Since the size of \textbf{c} is smaller than the size of \textbf{u}, this equation is underdetermined. To solve this problem, a Galerkin projection is used on the system:

ΦT(i=1NAfiA(μ)Ai)Φc(μ)=ΦTi=1Nbfib(μ)bi.\boldsymbol{\Phi}^T\left(\sum \limits_{i=1}^{N_A} f^A_i(\boldsymbol{\mu})\textbf{A}_i\right)\boldsymbol{\Phi}\textbf{c}(\boldsymbol{\mu}) = \boldsymbol{\Phi}^T\sum \limits_{i=1}^{N_b} f^b_i(\boldsymbol{\mu})\textbf{b}_i.(9)

By pulling the basis matrices into the summation, the following form is obtained.

(i=1NAfiA(μ)ΦTAiΦ)c(μ)=i=1Nbfib(μ)ΦTbi,\left(\sum \limits_{i=1}^{N_A} f^A_i(\boldsymbol{\mu})\boldsymbol{\Phi}^T\textbf{A}_i\boldsymbol{\Phi}\right)\textbf{c}(\boldsymbol{\mu}) = \sum \limits_{i=1}^{N_b} f^b_i(\boldsymbol{\mu})\boldsymbol{\Phi}^T\textbf{b}_i,(10)

where the reduced operators Air=ΦTAiΦ\textbf{A}_i^r=\boldsymbol{\Phi}^T\textbf{A}_i\boldsymbol{\Phi} and source terms bir=ΦTbi\textbf{b}_i^r=\boldsymbol{\Phi}^T\textbf{b}_i can be precomputed once the basis functions are available. Therefore, the reduced equation system that is solved to obtain c(μ)\textbf{c}(\boldsymbol{\mu}) is

(i=1NAfiA(μ)Air)c(μ)=i=1Nbfib(μ)bir.\left(\sum \limits_{i=1}^{N_A} f^A_i(\boldsymbol{\mu})\textbf{A}_i^r\right)\textbf{c}(\boldsymbol{\mu}) = \sum \limits_{i=1}^{N_b} f^b_i(\boldsymbol{\mu})\textbf{b}_i^r.(11)

It is important to note that this equation system is of size i=1Nvri×i=1Nvri\sum\limits_{i=1}^{N_v}r_i \times \sum\limits_{i=1}^{N_v}r_i, therefore it can be solved much faster than the original full-order system which is of size i=1NvNi×i=1NvNi\sum\limits_{i=1}^{N_v}N_i \times \sum\limits_{i=1}^{N_v}N_i. The trainer object only assembles the reduced operators and source terms, which are then transferred to PODReducedBasisSurrogate that is responsible for assembling and solving Eq. (11) and reconstructing the approximate solutions using Eq. (7) if needed.

The computation of the reduced operators consists of two step in the current implementation:

  1. Computing the effect of the full-order operator on the global basis functions: this step includes the creation of AiΦ\textbf{A}_i\boldsymbol{\Phi}. In practice, this is done by plugging in the basis function into a PODFullSolveMultiApp object which evaluates the residual for a given vector tag (defined using the "tag_names" input argument). The tagged residual is then transferred back to the trainer using a PODResidualTransfer object. In case when the residual from a kernel contains contributions to both the system matrix and the source term (e.g. Dirichlet BC or time derivative), certain input-level tricks can be used to separate these.

  2. Projection of the residual vectors: this step consists of computing the ΦT(AiΦ)\boldsymbol{\Phi}^T(\textbf{A}_i\boldsymbol{\Phi}) inner products.

As a final note, it must emphasized that even though obtaining snapshots and creating reduced operators is a computationally expensive procedure, it has to be carried out only once. After this initial investment every new evaluation for a new parameter sample involves the summation and scaling of small dense matrices, which is of low computational cost.

Note on parallelism

As of now, the training phase is implemented in a semi-parallel manner. This means that the snapshot generation, correlation matrix generation, base generation and the computation of the reduced operators are all executed in parallel. However, the eigenvalues and eigenvectors of the correlation matrices are obtained in serial. Therefore, this phase may experience considerable slowdown when the number of snapshots is large (above ~2000).

Example Input File Syntax

To get the snapshots, four essential blocks have to be added the main input file. First, a sampler has to be defined in Samplers to generate realizations for μ\boldsymbol{\mu}. These samples are then fed into a PODFullSolveMultiApp (defined in the MultiApps block) which is capable of running simulations for each parameter. The solution vectors (snapshots) for each run are added to the trainer by a PODSamplerSolutionTransfer defined in the Transfers block. It is important to mention that the multiapp and transfer objects need to know about the trainer. This can be ensured using the "trainer_name" input argument. The number of collected snapshots is defined in the sampler object using the "num_rows" parameter.

[Samplers]
  [sample]
    type = LatinHypercube
    distributions = 'k_dist alpha_dist S_dist'
    num_rows = 3
    execute_on = PRE_MULTIAPP_SETUP
    max_procs_per_row = 1
  []
[]
(contrib/moose/modules/stochastic_tools/test/tests/surrogates/pod_rb/internal/trainer.i)
[MultiApps]
  [sub]
    type = PODFullSolveMultiApp
    input_files = sub.i
    sampler = sample
    trainer_name = 'pod_rb'
    execute_on = 'timestep_begin final'
    max_procs_per_app = 1
  []
[]
(contrib/moose/modules/stochastic_tools/test/tests/surrogates/pod_rb/internal/trainer.i)

The global basis functions are then plugged back into the same PODFullSolveMultiApp by a new PODSamplerSolutionTransfer and the residuals are evaluated. The residuals are transferred back to the trainer using a PODResidualTransfer.

[Transfers]
  [param]
    type = SamplerParameterTransfer
    to_multi_app = sub
    sampler = sample
    parameters = 'Materials/k/prop_values Materials/alpha/prop_values Kernels/source/value'
    execute_on = 'timestep_begin'
    check_multiapp_execute_on = false
  []
  [snapshots]
    type = PODSamplerSolutionTransfer
    from_multi_app = sub
    sampler = sample
    trainer_name = 'pod_rb'
    execute_on = 'timestep_begin'
    check_multiapp_execute_on = false
  []
  [pod_modes]
    type = PODSamplerSolutionTransfer
    to_multi_app = sub
    sampler = sample
    trainer_name = 'pod_rb'
    execute_on = 'final'
    check_multiapp_execute_on = false
  []
  [res]
    type = PODResidualTransfer
    from_multi_app = sub
    sampler = sample
    trainer_name = "pod_rb"
    execute_on = 'final'
    check_multiapp_execute_on = false
  []
[]
(contrib/moose/modules/stochastic_tools/test/tests/surrogates/pod_rb/internal/trainer.i)

Finally, the PODReducedBasisTrainer is defined in the Trainers block. It requires the names of the variables ("var_names") one wishes to create reduced basis for and the names of the tags ("tag_names") associated with the affine components of the full-order operator. Additionally, a vector specifying which tag corresponds to a linear operator and which to a source term needs to be added as well ("tag_types"). The available tag types are:

Tag typeDescription
opOperator that acts on the solution vector (i.e. a matrix).
srcOperator that does not act on the solution vector (i.e. a source vector).
op_dirDirichlet operator that acts on the solution vector (i.e. Dirichlet penalty matrix).
src_dirDirichlet operator that does not act on the solution vector (i.e. Dirichlet penalty source vector).

To specify the allowed error in the snapshot reconstruction for reach variable (τi\tau_i), one needs to specify "error_res". In the following example the "execute_on" parameter is set to 'timestep_begin final' which means that at the beginning of the training it will create the basis functions while at the end it will create the reduced operators.

[Trainers]
  [pod_rb]
    type = PODReducedBasisTrainer
    var_names = 'u'
    error_res = '1e-9'
    tag_names = 'diff react bodyf'
    tag_types = 'op op src'
    execute_on = 'timestep_begin final'
  []
[]
(contrib/moose/modules/stochastic_tools/test/tests/surrogates/pod_rb/internal/trainer.i)

Input Parameters

  • error_resThe errors allowed in the snapshot reconstruction.

    C++ Type:std::vector<double>

    Unit:(no unit assumed)

    Controllable:No

    Description:The errors allowed in the snapshot reconstruction.

  • tag_namesNames of tags for the reduced operators.

    C++ Type:std::vector<std::string>

    Unit:(no unit assumed)

    Controllable:No

    Description:Names of tags for the reduced operators.

  • tag_typesList of keywords describing if the tags correspond to independent operators or not. (op/op_dir/src/src_dir)

    C++ Type:std::vector<std::string>

    Unit:(no unit assumed)

    Controllable:No

    Description:List of keywords describing if the tags correspond to independent operators or not. (op/op_dir/src/src_dir)

  • var_namesNames of variables we want to extract from solution vectors.

    C++ Type:std::vector<std::string>

    Unit:(no unit assumed)

    Controllable:No

    Description:Names of variables we want to extract from solution vectors.

Required Parameters

  • execute_onTIMESTEP_ENDThe list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html.

    Default:TIMESTEP_END

    C++ Type:ExecFlagEnum

    Unit:(no unit assumed)

    Options:NONE, INITIAL, LINEAR, NONLINEAR_CONVERGENCE, NONLINEAR, POSTCHECK, TIMESTEP_END, TIMESTEP_BEGIN, MULTIAPP_FIXED_POINT_END, MULTIAPP_FIXED_POINT_BEGIN, FINAL, CUSTOM

    Controllable:No

    Description:The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html.

  • filenameThe name of the file which will be associated with the saved/loaded data.

    C++ Type:FileName

    Unit:(no unit assumed)

    Controllable:No

    Description:The name of the file which will be associated with the saved/loaded data.

  • filenamesFiles where the eigenvalues are printed for each variable (if given).

    C++ Type:std::vector<std::string>

    Unit:(no unit assumed)

    Controllable:No

    Description:Files where the eigenvalues are printed for each variable (if given).

  • prop_getter_suffixAn optional suffix parameter that can be appended to any attempt to retrieve/get material properties. The suffix will be prepended with a '_' character.

    C++ Type:MaterialPropertyName

    Unit:(no unit assumed)

    Controllable:No

    Description:An optional suffix parameter that can be appended to any attempt to retrieve/get material properties. The suffix will be prepended with a '_' character.

  • use_interpolated_stateFalseFor the old and older state use projected material properties interpolated at the quadrature points. To set up projection use the ProjectedStatefulMaterialStorageAction.

    Default:False

    C++ Type:bool

    Unit:(no unit assumed)

    Controllable:No

    Description:For the old and older state use projected material properties interpolated at the quadrature points. To set up projection use the ProjectedStatefulMaterialStorageAction.

Optional Parameters

  • allow_duplicate_execution_on_initialFalseIn the case where this UserObject is depended upon by an initial condition, allow it to be executed twice during the initial setup (once before the IC and again after mesh adaptivity (if applicable).

    Default:False

    C++ Type:bool

    Unit:(no unit assumed)

    Controllable:No

    Description:In the case where this UserObject is depended upon by an initial condition, allow it to be executed twice during the initial setup (once before the IC and again after mesh adaptivity (if applicable).

  • control_tagsAdds user-defined labels for accessing object parameters via control logic.

    C++ Type:std::vector<std::string>

    Unit:(no unit assumed)

    Controllable:No

    Description:Adds user-defined labels for accessing object parameters via control logic.

  • enableTrueSet the enabled status of the MooseObject.

    Default:True

    C++ Type:bool

    Unit:(no unit assumed)

    Controllable:Yes

    Description:Set the enabled status of the MooseObject.

  • execution_order_group0Execution order groups are executed in increasing order (e.g., the lowest number is executed first). Note that negative group numbers may be used to execute groups before the default (0) group. Please refer to the user object documentation for ordering of user object execution within a group.

    Default:0

    C++ Type:int

    Unit:(no unit assumed)

    Controllable:No

    Description:Execution order groups are executed in increasing order (e.g., the lowest number is executed first). Note that negative group numbers may be used to execute groups before the default (0) group. Please refer to the user object documentation for ordering of user object execution within a group.

  • force_postauxFalseForces the UserObject to be executed in POSTAUX

    Default:False

    C++ Type:bool

    Unit:(no unit assumed)

    Controllable:No

    Description:Forces the UserObject to be executed in POSTAUX

  • force_preauxFalseForces the UserObject to be executed in PREAUX

    Default:False

    C++ Type:bool

    Unit:(no unit assumed)

    Controllable:No

    Description:Forces the UserObject to be executed in PREAUX

  • force_preicFalseForces the UserObject to be executed in PREIC during initial setup

    Default:False

    C++ Type:bool

    Unit:(no unit assumed)

    Controllable:No

    Description:Forces the UserObject to be executed in PREIC during initial setup

  • use_displaced_meshFalseWhether or not this object should use the displaced mesh for computation. Note that in the case this is true but no displacements are provided in the Mesh block the undisplaced mesh will still be used.

    Default:False

    C++ Type:bool

    Unit:(no unit assumed)

    Controllable:No

    Description:Whether or not this object should use the displaced mesh for computation. Note that in the case this is true but no displacements are provided in the Mesh block the undisplaced mesh will still be used.

Advanced Parameters

References

  1. René Pinnau. Model Reduction via Proper Orthogonal Decomposition, chapter 5, pages 95–109. Springer Berlin Heidelberg, Berlin, Heidelberg, 2008. doi:10.1007/978-3-540-78841-6_5.[BibTeX]