Using Stochastic Tools for the full-field reconstruction of multiphysics problems

The purpose of this example is to present a multiphysics model order reduction case using the Stochastic Tools Module. The problem of interest is a melt pool model using the Navier-Stokes module. The problem involves multiple uncertain process parameters and are interested in predicting full field quantities for unseen process parameter combinations using non-intrusive surrogate models in combination with Proper Orthogonal Decomposition (POD).

Problem Description

The problem of interest involves a transient 2D melt pool simulation for a 316L stainless steel with a moving Gaussian laser beam. The laser source melts the top layer of the metal and based on the equlibrium of the forces (surface tension, vapor recoil pressure) the surface of the melt pool can deform considerably. The changing geometry is taken into account using an Arbitrary Legrangian Eulerian (ALE) approach. The setup of the problem is depicted in Figure 1.

Figure 1: Problem setup for the laser melt pool simulation

The temperature-dependent thermal properties of the problem have been taken from the references listed in Table 1.

Table 1: Material properties for 316L Stailess steel

Property	Source
Dynamic viscosity	Kim (1975) $var element = document.getElementById("moose-equation-fc01e33e-3dce-41bd-9d20-0f29263a7a0a");katex.render("^\\dagger", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$
Density	Kim (1975)
Specific heat	Kim (1975)
Thermal conductivity	Pichler et al. (2022)
Surface tension	Pichler et al. (2022)
Vapor recoil pressure	Chen et al. (2021)

$var element = document.getElementById("moose-equation-2088caba-1b69-4c46-be01-91146f7ae4c1");katex.render("\\dagger", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ In the treatment of dynamic viscosity, following the recommendations in Noble et al. (2007), we artificially cap the minimum value to stabilize the fluid flow solver.

The input file used for the simulations is presented below:

Listing 1: Laser meltpool model input file

!include parameters.i

!include mesh.i

[Variables]
  [vel]
    family = LAGRANGE_VEC
  []
  [T]
  []
  [p]
  []
  [disp_x]
  []
  [disp_y]
  []
[]

!include physics_objects.i

[Executioner]
  type = Transient
  end_time = ${endtime}
  dtmin = 1e-10
  dtmax = 1e-5
  petsc_options_iname = '-pc_type -pc_factor_shift_type'
  petsc_options_value = 'lu       NONZERO'
  petsc_options = '-snes_converged_reason -ksp_converged_reason -options_left'
  solve_type = 'NEWTON'
  line_search = 'none'
  nl_max_its = 5
  l_max_its = 100
  [TimeStepper]
    type = IterationAdaptiveDT
    optimal_iterations = 5
    iteration_window = 1
    dt = ${timestep}
    linear_iteration_ratio = 1e6
    growth_factor = 1.25
  []
[]

[Reporters]
  [solution_storage]
    type = SolutionContainer
    execute_on = 'FINAL'
  []
[]

[Debug]
  show_var_residual_norms = true
[]

We see that it utilizes the !include syntax to read the process parameters togheter with the mesh and kernel and boundary condition objects from other input files. This input file design reduces duplication when it comes to the testing of the reduced-order models in subsequent sections.

Listing 2: Input file for specifying process parameters

# Process parameters
scanning_speed=1.0 # m/s
power=74 # W (this is the effective power so multiplied by eta)
R=1.25e-4 # m (this is the effective radius)

# Geometric parameters
thickness=1e-4 # m
xmin=-0.1e-3 # m
xmax=0.75e-3 # m
ymin=${fparse -thickness}
surfacetemp=300 # K (temperature at the other side of the plate)

# Time stepping parameters
endtime=4e-4 # s
timestep=${fparse endtime/1000} # s

Listing 3: Input file for specifying the mesh

[Mesh<<<{"href": "../../../syntax/Mesh/index.html"}>>>]
  [cmg]
    type = GeneratedMeshGenerator<<<{"description": "Create a line, square, or cube mesh with uniformly spaced or biased elements.", "href": "../../../source/meshgenerators/GeneratedMeshGenerator.html"}>>>
    dim<<<{"description": "The dimension of the mesh to be generated"}>>> = 2
    xmin<<<{"description": "Lower X Coordinate of the generated mesh"}>>> = ${xmin}
    xmax<<<{"description": "Upper X Coordinate of the generated mesh"}>>> = ${xmax}
    ymin<<<{"description": "Lower Y Coordinate of the generated mesh"}>>> = ${fparse ymin}
    ymax<<<{"description": "Upper Y Coordinate of the generated mesh"}>>> = 0
    nx<<<{"description": "Number of elements in the X direction"}>>> = 161
    ny<<<{"description": "Number of elements in the Y direction"}>>> = 50
  []
  displacements<<<{"description": "The variables corresponding to the x y z displacements of the mesh.  If this is provided then the displacements will be taken into account during the computation. Creation of the displaced mesh can be suppressed even if this is set by setting 'use_displaced_mesh = false'."}>>> = 'disp_x disp_y'
[]

Listing 4: Input file adding the objects that define the physics

[ICs<<<{"href": "../../../syntax/ICs/index.html"}>>>]
  [T]
    type = FunctionIC<<<{"description": "An initial condition that uses a normal function of x, y, z to produce values (and optionally gradients) for a field variable.", "href": "../../../source/ics/FunctionIC.html"}>>>
    variable<<<{"description": "The variable this initial condition is supposed to provide values for."}>>> = T
    function<<<{"description": "The initial condition function."}>>> = '(${surfacetemp} - 300) / ${thickness} * y + ${surfacetemp}'
  []
[]

[Kernels<<<{"href": "../../../syntax/Kernels/index.html"}>>>]
  [disp_x]
    type = Diffusion<<<{"description": "The Laplacian operator ($-\\nabla \\cdot \\nabla u$), with the weak form of $(\\nabla \\phi_i, \\nabla u_h)$.", "href": "../../../source/kernels/Diffusion.html"}>>>
    variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = disp_x
  []
  [disp_y]
    type = Diffusion<<<{"description": "The Laplacian operator ($-\\nabla \\cdot \\nabla u$), with the weak form of $(\\nabla \\phi_i, \\nabla u_h)$.", "href": "../../../source/kernels/Diffusion.html"}>>>
    variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = disp_y
  []
  [mass]
    type = INSADMass<<<{"description": "This class computes the mass equation residual and Jacobian contributions (the latter using automatic differentiation) for the incompressible Navier-Stokes equations.", "href": "../../../source/kernels/INSADMass.html"}>>>
    variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = p
    use_displaced_mesh<<<{"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."}>>> = true
  []
  [mass_pspg]
    type = INSADMassPSPG<<<{"description": "This class adds PSPG stabilization to the mass equation, enabling use of equal order shape functions for pressure and velocity variables", "href": "../../../source/kernels/INSADMassPSPG.html"}>>>
    variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = p
    use_displaced_mesh<<<{"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."}>>> = true
  []
  [momentum_time]
    type = INSADMomentumTimeDerivative<<<{"description": "This class computes the time derivative for the incompressible Navier-Stokes momentum equation.", "href": "../../../source/kernels/INSADMomentumTimeDerivative.html"}>>>
    variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = vel
    use_displaced_mesh<<<{"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."}>>> = true
  []
  [momentum_advection]
    type = INSADMomentumAdvection<<<{"description": "Adds the advective term to the INS momentum equation", "href": "../../../source/kernels/INSADMomentumAdvection.html"}>>>
    variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = vel
    use_displaced_mesh<<<{"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."}>>> = true
  []
  [momentum_mesh_advection]
    type = INSADMomentumMeshAdvection<<<{"description": "Corrects the convective derivative for situations in which the fluid mesh is dynamic.", "href": "../../../source/kernels/INSADMomentumMeshAdvection.html"}>>>
    variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = vel
    disp_x<<<{"description": "The x displacement"}>>> = disp_x
    disp_y<<<{"description": "The y displacement"}>>> = disp_y
    use_displaced_mesh<<<{"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."}>>> = true
  []
  [momentum_viscous]
    type = INSADMomentumViscous<<<{"description": "Adds the viscous term to the INS momentum equation", "href": "../../../source/kernels/INSADMomentumViscous.html"}>>>
    variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = vel
    use_displaced_mesh<<<{"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."}>>> = true
  []
  [momentum_pressure]
    type = INSADMomentumPressure<<<{"description": "Adds the pressure term to the INS momentum equation", "href": "../../../source/kernels/INSADMomentumPressure.html"}>>>
    variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = vel
    pressure<<<{"description": "The pressure"}>>> = p
    integrate_p_by_parts<<<{"description": "Whether to integrate the pressure term by parts"}>>> = true
    use_displaced_mesh<<<{"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."}>>> = true
  []
  [momentum_supg]
    type = INSADMomentumSUPG<<<{"description": "Adds the supg stabilization to the INS momentum equation", "href": "../../../source/kernels/INSADMomentumSUPG.html"}>>>
    variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = vel
    material_velocity<<<{"description": "A material property describing the velocity"}>>> = relative_velocity
    use_displaced_mesh<<<{"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."}>>> = true
  []
  [temperature_time]
    type = INSADHeatConductionTimeDerivative<<<{"description": "AD Time derivative term $\\rho c_p \\frac{\\partial T}{\\partial t}$ of the heat equation for quasi-constant specific heat $c_p$ and the density $\\rho$.", "href": "../../../source/kernels/INSADHeatConductionTimeDerivative.html"}>>>
    variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = T
    use_displaced_mesh<<<{"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."}>>> = true
  []
  [temperature_advection]
    type = INSADEnergyAdvection<<<{"description": "This class computes the residual and Jacobian contributions for temperature advection for a divergence free velocity field.", "href": "../../../source/kernels/INSADEnergyAdvection.html"}>>>
    variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = T
    use_displaced_mesh<<<{"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."}>>> = true
  []
  [temperature_mesh_advection]
    type = INSADEnergyMeshAdvection<<<{"description": "This class computes the residual and Jacobian contributions for temperature advection from mesh velocity in an ALE simulation.", "href": "../../../source/kernels/INSADEnergyMeshAdvection.html"}>>>
    variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = T
    disp_x<<<{"description": "The x displacement"}>>> = disp_x
    disp_y<<<{"description": "The y displacement"}>>> = disp_y
    use_displaced_mesh<<<{"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."}>>> = true
  []
  [temperature_conduction]
    type = ADHeatConduction<<<{"description": "Same as `Diffusion` in terms of physics/residual, but the Jacobian is computed using forward automatic differentiation", "href": "../../../source/kernels/ADHeatConduction.html"}>>>
    variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = T
    thermal_conductivity<<<{"description": "the name of the thermal conductivity material property"}>>> = 'k'
    use_displaced_mesh<<<{"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."}>>> = true
  []
  [temperature_supg]
    type = INSADEnergySUPG<<<{"description": "Adds the supg stabilization to the INS temperature/energy equation", "href": "../../../source/kernels/INSADEnergySUPG.html"}>>>
    variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = T
    velocity<<<{"description": "The velocity variable."}>>> = vel
    use_displaced_mesh<<<{"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."}>>> = true
  []
[]

[BCs<<<{"href": "../../../syntax/BCs/index.html"}>>>]
  [x_no_disp]
    type = DirichletBC<<<{"description": "Imposes the essential boundary condition $u=g$, where $g$ is a constant, controllable value.", "href": "../../../source/bcs/DirichletBC.html"}>>>
    variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = disp_x
    boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = 'bottom'
    value<<<{"description": "Value of the BC"}>>> = 0
  []
  [y_no_disp]
    type = DirichletBC<<<{"description": "Imposes the essential boundary condition $u=g$, where $g$ is a constant, controllable value.", "href": "../../../source/bcs/DirichletBC.html"}>>>
    variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = disp_y
    boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = 'bottom'
    value<<<{"description": "Value of the BC"}>>> = 0
  []
  [no_slip]
    type = ADVectorFunctionDirichletBC<<<{"description": "Imposes the essential boundary condition $\\vec{u}=\\vec{g}$, where $\\vec{g}$ components are calculated with functions.", "href": "../../../source/bcs/ADVectorFunctionDirichletBC.html"}>>>
    variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = vel
    boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = 'bottom right left'
  []
  [T_cold]
    type = DirichletBC<<<{"description": "Imposes the essential boundary condition $u=g$, where $g$ is a constant, controllable value.", "href": "../../../source/bcs/DirichletBC.html"}>>>
    variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = T
    boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = 'bottom'
    value<<<{"description": "Value of the BC"}>>> = 300
  []
  [radiation_flux]
    type = FunctionRadiativeBC<<<{"description": "Boundary condition for radiative heat exchange where the emissivity function is supplied by a Function.", "href": "../../../source/bcs/FunctionRadiativeBC.html"}>>>
    variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = T
    boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = 'top'
    emissivity_function<<<{"description": "Function describing emissivity for radiative boundary condition"}>>> = '1'
    Tinfinity<<<{"description": "Temperature of the body in radiative heat transfer."}>>> = 300
    stefan_boltzmann_constant<<<{"description": "The Stefan-Boltzmann constant."}>>> = 5.67e-8
    use_displaced_mesh<<<{"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."}>>> = true
  []
  [weld_flux]
    type = GaussianEnergyFluxBC<<<{"description": "Describes an incoming heat flux beam with a Gaussian profile", "href": "../../../source/bcs/GaussianEnergyFluxBC.html"}>>>
    variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = T
    boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = 'top'
    P0<<<{"description": "The total power of the beam."}>>> = ${power}
    R<<<{"description": "The radius at which the beam intensity falls to $1/e^2$ of its axis value."}>>> = ${R}
    x_beam_coord<<<{"description": "The x coordinate of the center of the beam as a function of time. Note that we will pass the origin as the spatial argument to the function; any spatial dependence in the passed-in function will be ignored"}>>> = '${scanning_speed}*t'
    y_beam_coord<<<{"description": "The y coordinate of the center of the beam as a function of time. Note that we will pass the origin as the spatial argument to the function; any spatial dependence in the passed-in function will be ignored"}>>> = '0'
    use_displaced_mesh<<<{"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."}>>> = true
  []
  [vapor_recoil]
    type = INSADVaporRecoilPressureMomentumFluxBC<<<{"description": "Vapor recoil pressure momentum flux", "href": "../../../source/bcs/INSADVaporRecoilPressureMomentumFluxBC.html"}>>>
    variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = vel
    boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = 'top'
    use_displaced_mesh<<<{"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."}>>> = true
  []
  [surface_tension]
    type = INSADSurfaceTensionBC<<<{"description": "Surface tension stresses.", "href": "../../../source/bcs/INSADSurfaceTensionBC.html"}>>>
    variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = vel
    boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = 'top'
    use_displaced_mesh<<<{"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."}>>> = true
  []
  [displace_x_top]
    type = INSADDisplaceBoundaryBC<<<{"description": "Boundary condition for displacing a boundary", "href": "../../../source/bcs/INSADDisplaceBoundaryBC.html"}>>>
    boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = 'top'
    variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = 'disp_x'
    velocity<<<{"description": "The velocity at which to displace. A functor is any of the following: a variable, a functor material property, a function, a postprocessor or a number."}>>> = 'vel'
    component<<<{"description": "What component of velocity/displacement this object is acting on."}>>> = 0
    associated_subdomain<<<{"description": "The subdomain that the boundary nodeset is associated with. This will be passed to the coupled functor for unambigious evaluation (e.g. at the edge of the node-patch where we might run into the intersection of subdomains"}>>> = 0
  []
  [displace_y_top]
    type = INSADDisplaceBoundaryBC<<<{"description": "Boundary condition for displacing a boundary", "href": "../../../source/bcs/INSADDisplaceBoundaryBC.html"}>>>
    boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = 'top'
    variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = 'disp_y'
    velocity<<<{"description": "The velocity at which to displace. A functor is any of the following: a variable, a functor material property, a function, a postprocessor or a number."}>>> = 'vel'
    component<<<{"description": "What component of velocity/displacement this object is acting on."}>>> = 1
    associated_subdomain<<<{"description": "The subdomain that the boundary nodeset is associated with. This will be passed to the coupled functor for unambigious evaluation (e.g. at the edge of the node-patch where we might run into the intersection of subdomains"}>>> = 0
  []
  [displace_x_top_dummy]
    type = INSADDummyDisplaceBoundaryIntegratedBC<<<{"description": "This object adds Jacobian entries for the boundary displacement dependence on the velocity", "href": "../../../source/bcs/INSADDummyDisplaceBoundaryIntegratedBC.html"}>>>
    boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = 'top'
    variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = 'disp_x'
    velocity<<<{"description": "The velocity at which to displace. A functor is any of the following: a variable, a functor material property, a function, a postprocessor or a number."}>>> = 'vel'
    component<<<{"description": "What component of velocity/displacement this object is acting on."}>>> = 0
  []
  [displace_y_top_dummy]
    type = INSADDummyDisplaceBoundaryIntegratedBC<<<{"description": "This object adds Jacobian entries for the boundary displacement dependence on the velocity", "href": "../../../source/bcs/INSADDummyDisplaceBoundaryIntegratedBC.html"}>>>
    boundary<<<{"description": "The list of boundary IDs from the mesh where this object applies"}>>> = 'top'
    variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = 'disp_y'
    velocity<<<{"description": "The velocity at which to displace. A functor is any of the following: a variable, a functor material property, a function, a postprocessor or a number."}>>> = 'vel'
    component<<<{"description": "What component of velocity/displacement this object is acting on."}>>> = 1
  []
[]

[Materials<<<{"href": "../../../syntax/Materials/index.html"}>>>]
  [ins_mat]
    type = INSADStabilized3Eqn<<<{"description": "This is the material class used to compute the stabilization parameter tau for momentum and tau_energy for the energy equation.", "href": "../../../source/materials/INSADStabilized3Eqn.html"}>>>
    velocity<<<{"description": "The velocity"}>>> = vel
    pressure<<<{"description": "The pressure"}>>> = p
    temperature<<<{"description": "The temperature"}>>> = T
    use_displaced_mesh<<<{"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."}>>> = true
  []
  [steel]
    type = LaserWeld316LStainlessSteel
    temperature = T
    use_displaced_mesh = true
  []
  [steel_boundary]
    type = LaserWeld316LStainlessSteelBoundary
    boundary = 'top'
    temperature = T
    use_displaced_mesh = true
  []
  [const]
    type = GenericConstantMaterial<<<{"description": "Declares material properties based on names and values prescribed by input parameters.", "href": "../../../source/materials/GenericConstantMaterial.html"}>>>
    prop_names<<<{"description": "The names of the properties this material will have"}>>> = 'abs sb_constant'
    prop_values<<<{"description": "The values associated with the named properties"}>>> = '1 5.67e-8'
    use_displaced_mesh<<<{"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."}>>> = true
  []
[]

Process parameters

A total of two process parameters have been identified to be changable: the effective laser power and the effective laser radius. They can change in the intervals shown in Table 2 assuming a uniform distribution.

Table 2: The distributions used for the process parameters

Property	Min.	Max.
Effective laser power	60 W	74 W
Effective laser radius	125 $var element = document.getElementById("moose-equation-3086b703-4644-4b2e-a59b-986787c62ffc");katex.render("\\mu m", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$	155 $var element = document.getElementById("moose-equation-44c0cd4c-391c-4839-8892-a1bc30a9c9fa");katex.render("\\mu m", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$

Quantities of Interest

The quantity of interest in our case is the whole temperature field at the end of the simulation (0.4~ms). To be able to reconstruct the whole field we utilize Proper Orthogonal Decomposition (POD) in conjunction with a Multi-output Gaussian Process (MOGP).

Results

Reference results

Two reference results are presented here: 1. Results for a simulation with low effective laser power ( $var element = document.getElementById("moose-equation-3dc3d7d4-b9d1-48bd-8f4a-2a9404b11421");katex.render("60~W", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ ) and wide effective laser radius ( $var element = document.getElementById("moose-equation-b68fad28-f02d-4f38-bf7a-81ac47fabd8c");katex.render("155~\\mu m", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ ) 2. Results for a simulation with high effective laser power ( $var element = document.getElementById("moose-equation-a2012897-d59c-472b-85ca-eb18d462db80");katex.render("74~W", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ ) and narrow effective laser radius ( $var element = document.getElementById("moose-equation-f6cfb74c-b983-4374-9a7f-4110efd94b08");katex.render("125~\\mu m", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ )

The reference results are filtered to show the melt pool shape. We see that for the first case, the steel melts, but no displacement is visible due to little to no evaporation. On the other hand, the increase in laser power results in an increased evaporation resulting in significant surface deformation.

Figure 2: Simulation results at the final time step with laser power ( $var element = document.getElementById("moose-equation-44d9ef7d-0abc-4689-97cb-d69e5909572b");katex.render("60~W", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ ) and wide effective laser radius ( $var element = document.getElementById("moose-equation-fedbd0d2-0970-49e2-9e80-666590d4e14d");katex.render("155~\\mu m", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ )

Figure 3: Simulation results at the final time step with laser power ( $var element = document.getElementById("moose-equation-51a09f23-365d-4491-a79a-9fd2d7c35039");katex.render("74~W", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ ) and narrow effective laser radius ( $var element = document.getElementById("moose-equation-523727f9-6947-4dfc-9217-3c0330d1daa2");katex.render("125~\\mu m", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ )

Training the reduced-order model

The reduced-order model was trained by first running the reference model with 45 different process parameter combinations in a Monte Carlo sampler. The input file used for training is presented below.

Listing 5:

[StochasticTools<<<{"href": "../../../syntax/StochasticTools/index.html"}>>>]
[]

[Distributions<<<{"href": "../../../syntax/Distributions/index.html"}>>>]
  [R_dist]
    type = Uniform<<<{"description": "Continuous uniform distribution.", "href": "../../../source/distributions/Uniform.html"}>>>
    lower_bound<<<{"description": "Distribution lower bound"}>>> = 1.25E-4
    upper_bound<<<{"description": "Distribution upper bound"}>>> = 1.55E-4
  []
  [power_dist]
    type = Uniform<<<{"description": "Continuous uniform distribution.", "href": "../../../source/distributions/Uniform.html"}>>>
    lower_bound<<<{"description": "Distribution lower bound"}>>> = 60
    upper_bound<<<{"description": "Distribution upper bound"}>>> = 74
  []
[]

[Samplers<<<{"href": "../../../syntax/Samplers/index.html"}>>>]
  [train]
    type = MonteCarlo<<<{"description": "Monte Carlo Sampler.", "href": "../../../source/samplers/MonteCarloSampler.html"}>>>
    num_rows<<<{"description": "The number of rows per matrix to generate."}>>> = 45
    distributions<<<{"description": "The distribution names to be sampled, the number of distributions provided defines the number of columns per matrix."}>>> = 'R_dist power_dist'
    execute_on<<<{"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."}>>> = PRE_MULTIAPP_SETUP
    min_procs_per_row<<<{"description": "This will ensure that the sampler is partitioned properly when 'MultiApp/*/min_procs_per_app' is specified. It is not recommended to use otherwise."}>>> = 2
    max_procs_per_row<<<{"description": "This will ensure that the sampler is partitioned properly when 'MultiApp/*/max_procs_per_app' is specified. It is not recommended to use otherwise."}>>> = 2
  []
[]

[MultiApps<<<{"href": "../../../syntax/MultiApps/index.html"}>>>]
  [worker]
    type = SamplerFullSolveMultiApp<<<{"description": "Creates a full-solve type sub-application for each row of each Sampler matrix.", "href": "../../../source/multiapps/SamplerFullSolveMultiApp.html"}>>>
    input_files<<<{"description": "The input file for each App.  If this parameter only contains one input file it will be used for all of the Apps.  When using 'positions_from_file' it is also admissable to provide one input_file per file."}>>> = 2d.i
    sampler<<<{"description": "The Sampler object to utilize for creating MultiApps."}>>> = train
    mode<<<{"description": "The operation mode, 'normal' creates one sub-application for each row in the Sampler and 'batch-reset' and 'batch-restore' creates N sub-applications, where N is the minimum of 'num_rows' in the Sampler and floor(number of processes / min_procs_per_app). To run the rows in the Sampler, 'batch-reset' will destroy and re-create sub-apps as needed, whereas the 'batch-restore' will backup and restore sub-apps to the initial state prior to execution, without destruction."}>>> = batch-reset
    min_procs_per_app<<<{"description": "Minimum number of processors to give to each App in this MultiApp.  Useful for larger, distributed mesh solves."}>>> = 2
    max_procs_per_app<<<{"description": "Maximum number of processors to give to each App in this MultiApp.  Useful for restricting small solves to just a few procs so they don't get spread out"}>>> = 2
  []
[]

[Controls<<<{"href": "../../../syntax/Controls/index.html"}>>>]
  [cmdline]
    type = MultiAppSamplerControl<<<{"description": "Control for modifying the command line arguments of MultiApps.", "href": "../../../source/controls/MultiAppSamplerControl.html"}>>>
    multi_app<<<{"description": "The name of the MultiApp to control."}>>> = worker
    sampler<<<{"description": "The Sampler object to utilize for altering the command line options of the MultiApp."}>>> = train
    param_names<<<{"description": "The names of the command line parameters to set via the sampled data."}>>> = 'R power'
  []
[]

[Transfers<<<{"href": "../../../syntax/Transfers/index.html"}>>>]
  [solution_transfer]
    type = SerializedSolutionTransfer<<<{"description": "Serializes and transfers solution vectors for given variables from sub-applications.", "href": "../../../source/transfers/SerializedSolutionTransfer.html"}>>>
    parallel_storage<<<{"description": "The name of the parallel storage reporter."}>>> = parallel_storage
    from_multi_app<<<{"description": "The name of the MultiApp to receive data from"}>>> = worker
    sampler<<<{"description": "A the Sampler object that Transfer is associated.."}>>> = train
    solution_container<<<{"description": "The name of the solution container on the subapp."}>>> = solution_storage
    variables<<<{"description": "The names of the variables which should be serialized and transferred to this application."}>>> = "T"
    serialize_on_root<<<{"description": "If we want to gather the solution fields only on the root processors of the subapps before transfering to the main app."}>>> = true
  []
[]

[Reporters<<<{"href": "../../../syntax/Reporters/index.html"}>>>]
  [parallel_storage]
    type = ParallelSolutionStorage<<<{"description": "Parallel container to store serialized solution fields from simulations on sub-applications.", "href": "../../../source/reporters/ParallelSolutionStorage.html"}>>>
    variables<<<{"description": "The names of the variables whose serialized solution this object is supposed to receive."}>>> = "T"
    outputs<<<{"description": "Vector of output names where you would like to restrict the output of variables(s) associated with this object"}>>> = none
  []
  [reduced_sol]
    type = MappingReporter<<<{"description": "A reporter which can map full solution fields to a latent space for given variables.", "href": "../../../source/reporters/MappingReporter.html"}>>>
    sampler<<<{"description": "Sampler be able to identify how the samples are distributed among the processes. Only needed if parallel storage is defined. It is important to have the same sampler here as the one used to prepare the snapshots in the parallel storage."}>>> = train
    parallel_storage<<<{"description": "The storage space where the snapshots are stored. These snapshots are used to build the mapping. If this parameter is not specified, the reporter will fetch the variable from the nonlinear system."}>>> = parallel_storage
    mapping<<<{"description": "Name of the mapping object."}>>> = pod
    variables<<<{"description": "The names of the variables which need to be mapped to the latent space."}>>> = "T"
    outputs<<<{"description": "Vector of output names where you would like to restrict the output of variables(s) associated with this object"}>>> = json
    execute_on<<<{"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."}>>> = final
    parallel_type<<<{"description": "This parameter will determine how the stochastic data is gathered. It is common for outputting purposes that this parameter be set to ROOT, otherwise, many files will be produced showing the values on each processor. However, if there are lot of samples, gathering on root may be memory restrictive."}>>> = ROOT
  []
  [matrix]
    type = StochasticMatrix<<<{"description": "Tool for extracting Sampler object data and storing data from stochastic simulations.", "href": "../../../source/reporters/StochasticMatrix.html"}>>>
    sampler<<<{"description": "The sample from which to extract distribution data."}>>> = train
    parallel_type<<<{"description": "This parameter will determine how the stochastic data is gathered. It is common for outputting purposes that this parameter be set to ROOT, otherwise, many files will be produced showing the values on each processor. However, if there are lot of samples, gathering on root may be memory restrictive."}>>> = ROOT
  []
  [svd]
    type = SingularTripletReporter<<<{"description": "Tool for accessing and outputting the singular triplets of a singular value decomposition in PODMapping.", "href": "../../../source/reporters/SingularTripletReporter.html"}>>>
    pod_mapping<<<{"description": "The PODMapping object whose singular triplets should be printed."}>>> = pod
    variables<<<{"description": "The names of the variables whose SVD should be printed."}>>> = "T"
    execute_on<<<{"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."}>>> = final
  []
[]

[VariableMappings<<<{"href": "../../../syntax/VariableMappings/index.html"}>>>]
  [pod]
    type = PODMapping<<<{"description": "Class which provides a Proper Orthogonal Decomposition-based mapping between full-order and reduced-order spaces.", "href": "../../../source/variablemappings/PODMapping.html"}>>>
    solution_storage<<<{"description": "The name of the storage reporter where the snapshots are located."}>>> = parallel_storage
    variables<<<{"description": "The names of the variables which need a mapping."}>>> = "T"
    num_modes_to_compute<<<{"description": "The number of modes that this object should compute. Modes with 0 eigenvalues are filtered out, so the real number of modes might be lower than this. This is also used for setting the subspace sizes for distributed singular value solves. By default, the subspace used for the SVD is twice as big as the number of requested vectors. For more information see the SLEPc manual. If not specified, only one mode is computed per variable."}>>> = '8'
    extra_slepc_options<<<{"description": "Additional options for the singular/eigenvalue solvers in SLEPc."}>>> = "-svd_monitor_all"
  []
[]

[Trainers<<<{"href": "../../../syntax/Trainers/index.html"}>>>]
  [mogp]
    type = GaussianProcessTrainer<<<{"description": "Provides data preperation and training for a single- or multi-output Gaussian Process surrogate model.", "href": "../../../source/trainers/GaussianProcessTrainer.html"}>>>
    execute_on<<<{"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."}>>> = final
    covariance_function<<<{"description": "Name of covariance function."}>>> = 'lmc'
    standardize_params<<<{"description": "Standardize (center and scale) training parameters (x values)"}>>> = 'true'
    standardize_data<<<{"description": "Standardize (center and scale) training data (y values)"}>>> = 'true'
    sampler<<<{"description": "Sampler used to create predictor and response data."}>>> = train
    response_type<<<{"description": "Response data type."}>>> = vector_real
    response<<<{"description": "Reporter value of response results, can be vpp with <vpp_name>/<vector_name> or sampler column with 'sampler/col_<index>'."}>>> = reduced_sol/T_pod
    tune_parameters<<<{"description": "Select hyperparameters to be tuned"}>>> = 'lmc:acoeff_0 lmc:lambda_0 covar:signal_variance covar:length_factor'
    tuning_min<<<{"description": "Minimum allowable tuning value"}>>> = '1e-9 1e-9 1e-9 1e-9'
    tuning_max<<<{"description": "Maximum allowable tuning value"}>>> = '1e16 1e16 1e16  1e16'
    num_iters<<<{"description": "Tolerance value for Adam optimization"}>>> = 10000
    learning_rate<<<{"description": "The learning rate for Adam optimization"}>>> = 0.0005
    show_every_nth_iteration<<<{"description": "Switch to show Adam optimization loss values at every nth step. If 0, nothing is showed."}>>> = 10
  []
[]

[Covariance<<<{"href": "../../../syntax/Covariance/index.html"}>>>]
  [covar]
    type = SquaredExponentialCovariance<<<{"description": "Squared Exponential covariance function.", "href": "../../../source/covariances/SquaredExponentialCovariance.html"}>>>
    signal_variance<<<{"description": "Signal Variance ($\\sigma_f^2$) to use for kernel calculation."}>>> = 1.0
    noise_variance<<<{"description": "Noise Variance ($\\sigma_n^2$) to use for kernel calculation."}>>> = 0.0
    length_factor<<<{"description": "Length factors to use for Covariance Kernel"}>>> = '0.1 0.1'
  []
  [lmc]
    type = LMC<<<{"description": "Covariance function for multioutput Gaussian Processes based on the Linear Model of Coregionalization (LMC).", "href": "../../../source/covariances/LMC.html"}>>>
    covariance_functions<<<{"description": "Covariance functions that this covariance function depends on."}>>> = covar
    num_outputs<<<{"description": "The number of outputs expected for this covariance function."}>>> = 8
    num_latent_funcs<<<{"description": "The number of latent functions for the expansion of the outputs."}>>> = 1
  []
[]

[Outputs<<<{"href": "../../../syntax/Outputs/index.html"}>>>]
  [json]
    type = JSON<<<{"description": "Output for Reporter values using JSON format.", "href": "../../../source/outputs/JSONOutput.html"}>>>
    execute_on<<<{"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."}>>> = FINAL
    execute_system_information_on<<<{"description": "Control when the output of the simulation information occurs"}>>> = NONE
  []
  [pod_out]
    type = MappingOutput<<<{"description": "Output for mapping model data.", "href": "../../../source/outputs/MappingOutput.html"}>>>
    mappings<<<{"description": "A list of Mapping objects to output."}>>> = pod
    execute_on<<<{"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."}>>> = FINAL
  []
  [mogp_out]
    type = SurrogateTrainerOutput<<<{"description": "Output for trained surrogate model data.", "href": "../../../source/outputs/SurrogateTrainerOutput.html"}>>>
    trainers<<<{"description": "A list of SurrogateTrainer objects to output."}>>> = mogp
    execute_on<<<{"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."}>>> = FINAL
  []
[]

The parts that require special attention are the snapshot gathering parts in the transfer and reporter blocks.

Listing 6: Transfers and Reporters needed for snapshot collection

[Transfers<<<{"href": "../../../syntax/Transfers/index.html"}>>>]
  [solution_transfer]
    type = SerializedSolutionTransfer<<<{"description": "Serializes and transfers solution vectors for given variables from sub-applications.", "href": "../../../source/transfers/SerializedSolutionTransfer.html"}>>>
    parallel_storage<<<{"description": "The name of the parallel storage reporter."}>>> = parallel_storage
    from_multi_app<<<{"description": "The name of the MultiApp to receive data from"}>>> = worker
    sampler<<<{"description": "A the Sampler object that Transfer is associated.."}>>> = train
    solution_container<<<{"description": "The name of the solution container on the subapp."}>>> = solution_storage
    variables<<<{"description": "The names of the variables which should be serialized and transferred to this application."}>>> = "T"
    serialize_on_root<<<{"description": "If we want to gather the solution fields only on the root processors of the subapps before transfering to the main app."}>>> = true
  []
[]

[Reporters<<<{"href": "../../../syntax/Reporters/index.html"}>>>]
  [parallel_storage]
    type = ParallelSolutionStorage<<<{"description": "Parallel container to store serialized solution fields from simulations on sub-applications.", "href": "../../../source/reporters/ParallelSolutionStorage.html"}>>>
    variables<<<{"description": "The names of the variables whose serialized solution this object is supposed to receive."}>>> = "T"
    outputs<<<{"description": "Vector of output names where you would like to restrict the output of variables(s) associated with this object"}>>> = none
  []
  [reduced_sol]
    type = MappingReporter<<<{"description": "A reporter which can map full solution fields to a latent space for given variables.", "href": "../../../source/reporters/MappingReporter.html"}>>>
    sampler<<<{"description": "Sampler be able to identify how the samples are distributed among the processes. Only needed if parallel storage is defined. It is important to have the same sampler here as the one used to prepare the snapshots in the parallel storage."}>>> = train
    parallel_storage<<<{"description": "The storage space where the snapshots are stored. These snapshots are used to build the mapping. If this parameter is not specified, the reporter will fetch the variable from the nonlinear system."}>>> = parallel_storage
    mapping<<<{"description": "Name of the mapping object."}>>> = pod
    variables<<<{"description": "The names of the variables which need to be mapped to the latent space."}>>> = "T"
    outputs<<<{"description": "Vector of output names where you would like to restrict the output of variables(s) associated with this object"}>>> = json
    execute_on<<<{"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."}>>> = final
    parallel_type<<<{"description": "This parameter will determine how the stochastic data is gathered. It is common for outputting purposes that this parameter be set to ROOT, otherwise, many files will be produced showing the values on each processor. However, if there are lot of samples, gathering on root may be memory restrictive."}>>> = ROOT
  []
[]

[VariableMappings<<<{"href": "../../../syntax/VariableMappings/index.html"}>>>]
  [pod]
    type = PODMapping<<<{"description": "Class which provides a Proper Orthogonal Decomposition-based mapping between full-order and reduced-order spaces.", "href": "../../../source/variablemappings/PODMapping.html"}>>>
    solution_storage<<<{"description": "The name of the storage reporter where the snapshots are located."}>>> = parallel_storage
    variables<<<{"description": "The names of the variables which need a mapping."}>>> = "T"
    num_modes_to_compute<<<{"description": "The number of modes that this object should compute. Modes with 0 eigenvalues are filtered out, so the real number of modes might be lower than this. This is also used for setting the subspace sizes for distributed singular value solves. By default, the subspace used for the SVD is twice as big as the number of requested vectors. For more information see the SLEPc manual. If not specified, only one mode is computed per variable."}>>> = '8'
    extra_slepc_options<<<{"description": "Additional options for the singular/eigenvalue solvers in SLEPc."}>>> = "-svd_monitor_all"
  []
[]

With the snapshots all available in one place, we can extract the POD modes necessary for the dimensionality reduction. This is done using the PODMapping reporter which not only creates the decomposition but also maps the created snapshots onto a low-dimensional latent space.

Figure 4: Scree plot of the first 8 singular values extracted from the training data.

In this specific example, the dimensionality of the latent space was selected to be 8. Once the coordinates of the snapshots are available in the latent space, we can fit a MOGP on them. This is done in the Trainers block with the corresponding covariance functions defined in the Covariance block.

Listing 7: Definition of the Gaussian Process surrogate used for fitting the coordinates in the latent space.

[Trainers<<<{"href": "../../../syntax/Trainers/index.html"}>>>]
  [mogp]
    type = GaussianProcessTrainer<<<{"description": "Provides data preperation and training for a single- or multi-output Gaussian Process surrogate model.", "href": "../../../source/trainers/GaussianProcessTrainer.html"}>>>
    execute_on<<<{"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."}>>> = final
    covariance_function<<<{"description": "Name of covariance function."}>>> = 'lmc'
    standardize_params<<<{"description": "Standardize (center and scale) training parameters (x values)"}>>> = 'true'
    standardize_data<<<{"description": "Standardize (center and scale) training data (y values)"}>>> = 'true'
    sampler<<<{"description": "Sampler used to create predictor and response data."}>>> = train
    response_type<<<{"description": "Response data type."}>>> = vector_real
    response<<<{"description": "Reporter value of response results, can be vpp with <vpp_name>/<vector_name> or sampler column with 'sampler/col_<index>'."}>>> = reduced_sol/T_pod
    tune_parameters<<<{"description": "Select hyperparameters to be tuned"}>>> = 'lmc:acoeff_0 lmc:lambda_0 covar:signal_variance covar:length_factor'
    tuning_min<<<{"description": "Minimum allowable tuning value"}>>> = '1e-9 1e-9 1e-9 1e-9'
    tuning_max<<<{"description": "Maximum allowable tuning value"}>>> = '1e16 1e16 1e16  1e16'
    num_iters<<<{"description": "Tolerance value for Adam optimization"}>>> = 10000
    learning_rate<<<{"description": "The learning rate for Adam optimization"}>>> = 0.0005
    show_every_nth_iteration<<<{"description": "Switch to show Adam optimization loss values at every nth step. If 0, nothing is showed."}>>> = 10
  []
[]

[Covariance<<<{"href": "../../../syntax/Covariance/index.html"}>>>]
  [covar]
    type = SquaredExponentialCovariance<<<{"description": "Squared Exponential covariance function.", "href": "../../../source/covariances/SquaredExponentialCovariance.html"}>>>
    signal_variance<<<{"description": "Signal Variance ($\\sigma_f^2$) to use for kernel calculation."}>>> = 1.0
    noise_variance<<<{"description": "Noise Variance ($\\sigma_n^2$) to use for kernel calculation."}>>> = 0.0
    length_factor<<<{"description": "Length factors to use for Covariance Kernel"}>>> = '0.1 0.1'
  []
  [lmc]
    type = LMC<<<{"description": "Covariance function for multioutput Gaussian Processes based on the Linear Model of Coregionalization (LMC).", "href": "../../../source/covariances/LMC.html"}>>>
    covariance_functions<<<{"description": "Covariance functions that this covariance function depends on."}>>> = covar
    num_outputs<<<{"description": "The number of outputs expected for this covariance function."}>>> = 8
    num_latent_funcs<<<{"description": "The number of latent functions for the expansion of the outputs."}>>> = 1
  []
[]

One can then save both the POD mapping and MOGP into restratable structures for later usage.

[Outputs<<<{"href": "../../../syntax/Outputs/index.html"}>>>]
  [json]
    type = JSON<<<{"description": "Output for Reporter values using JSON format.", "href": "../../../source/outputs/JSONOutput.html"}>>>
    execute_on<<<{"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."}>>> = FINAL
    execute_system_information_on<<<{"description": "Control when the output of the simulation information occurs"}>>> = NONE
  []
  [pod_out]
    type = MappingOutput<<<{"description": "Output for mapping model data.", "href": "../../../source/outputs/MappingOutput.html"}>>>
    mappings<<<{"description": "A list of Mapping objects to output."}>>> = pod
    execute_on<<<{"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."}>>> = FINAL
  []
  [mogp_out]
    type = SurrogateTrainerOutput<<<{"description": "Output for trained surrogate model data.", "href": "../../../source/outputs/SurrogateTrainerOutput.html"}>>>
    trainers<<<{"description": "A list of SurrogateTrainer objects to output."}>>> = mogp
    execute_on<<<{"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."}>>> = FINAL
  []
[]

Evaluating the error of the ROM

We evaluate the ROM by comparing reconstructed temperature fields at unseen process parameter combinations to full-order model results at the same samples. The metric for comparison is the relative $var element = document.getElementById("moose-equation-377bc434-7c49-4e2e-9f7b-9bf41f157ff3");katex.render("L^2", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ error defined as:

$var element = document.getElementById("moose-equation-7b8a15bd-74b3-4712-a381-efbee2c4c280");katex.render("e = \\frac{||T_{ROM}-T_{FOM}||}{||T_{FOM}||}", element, {displayMode:true,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$

This error is generated in a new input file, which computes the reconstructed field and solves the full-order model at the same time.

Listing 8: Test input file for comparing the reconstructed fields to the reference solutions

!include parameters.i

!include mesh.i

[Variables]
  [vel]
    family = LAGRANGE_VEC
  []
  [T]
  []
  [p]
  []
  [disp_x]
  []
  [disp_y]
  []
[]

[AuxVariables]
  [T_reconst]
  []
[]

!include physics_objects.i

[UserObjects]
  [inverse]
    type = InverseMapping
    mapping = pod
    variable_to_fill = "T_reconst"
    variable_to_reconstruct = "T"
    surrogate = mogp
    parameters = '${R} ${power}'
    execute_on = TIMESTEP_END
  []
[]

[Surrogates]
  [mogp]
    type = GaussianProcessSurrogate
    filename = "train_mogp_out_mogp.rd"
  []
[]

[VariableMappings]
  [pod]
    type = PODMapping
    filename = "train_pod_out_pod.rd"
  []
[]

[Executioner]
  type = Transient
  end_time = ${endtime}
  dtmin = 1e-10
  dtmax = 1e-5
  petsc_options_iname = '-pc_type -pc_factor_shift_type'
  petsc_options_value = 'lu       NONZERO'
  petsc_options = '-snes_converged_reason -ksp_converged_reason -options_left'
  solve_type = 'NEWTON'
  line_search = 'none'
  nl_max_its = 16
  l_max_its = 100
  [TimeStepper]
    type = IterationAdaptiveDT
    optimal_iterations = 10
    iteration_window = 4
    dt = ${timestep}
    linear_iteration_ratio = 1e6
    growth_factor = 1.25
  []
[]

[Postprocessors]
  [l2error]
    type = ElementL2Difference
    variable = T
    other_variable = T_reconst
  []
[]

[Debug]
  show_var_residual_norms = true
[]

[Postprocessors]
[]

We see that duplication is minimized by reusing the same building blocks as the original input file. In this case, we load both the MOGP and the POD mapping from a file:

[VariableMappings<<<{"href": "../../../syntax/VariableMappings/index.html"}>>>]
  [pod]
    type = PODMapping<<<{"description": "Class which provides a Proper Orthogonal Decomposition-based mapping between full-order and reduced-order spaces.", "href": "../../../source/variablemappings/PODMapping.html"}>>>
    filename<<<{"description": "The name of the file which will be associated with the saved/loaded data."}>>> = "train_pod_out_pod.rd"
  []
[]

[Surrogates<<<{"href": "../../../syntax/Surrogates/index.html"}>>>]
  [mogp]
    type = GaussianProcessSurrogate<<<{"description": "Computes and evaluates Gaussian Process surrogate model.", "href": "../../../source/surrogates/GaussianProcessSurrogate.html"}>>>
    filename<<<{"description": "The name of the file which will be associated with the saved/loaded data."}>>> = "train_mogp_out_mogp.rd"
  []
[]

And use them in an object that reconstructs the field variable and inserts it into a T_reconst auxiliary variable.

[UserObjects<<<{"href": "../../../syntax/UserObjects/index.html"}>>>]
  [inverse]
    type = InverseMapping<<<{"description": "Evaluates surrogate models and maps the results back to a full solution field for given variables.", "href": "../../../source/userobjects/InverseMapping.html"}>>>
    mapping<<<{"description": "The name of the mapping object which provides the inverse mapping function."}>>> = pod
    variable_to_fill<<<{"description": "The names of the variables that this object is supposed to populate with the reconstructed results."}>>> = "T_reconst"
    variable_to_reconstruct<<<{"description": "The names of the variables in the nonlinear system which we would like to approximate. This is important for DoF information."}>>> = "T"
    surrogate<<<{"description": "The names of the surrogates for each variable."}>>> = mogp
    parameters<<<{"description": "The input parameters for the surrogate. If no surrogate is supplied these are assumed to be the coordinates in the latent space."}>>> = '${R} ${power}'
    execute_on<<<{"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."}>>> = TIMESTEP_END
  []
[]

Then, the integrated L2 error is computed by the appropriate postprocessor.

[Postprocessors<<<{"href": "../../../syntax/Postprocessors/index.html"}>>>]
  [l2error]
    type = ElementL2Difference<<<{"description": "Computes the element-wise L2 difference between the current variable and a coupled variable.", "href": "../../../source/postprocessors/ElementL2Difference.html"}>>>
    variable<<<{"description": "The name of the variable that this object operates on"}>>> = T
    other_variable<<<{"description": "The variable to compare to"}>>> = T_reconst
  []
[]

This input file is run for 90 samples of new process parameter combination selected by a Monte Carlo sampler. The input file for automating the process of error analysis is presented below.

Listing 9: Input file for the automation of the testing phase.

[StochasticTools<<<{"href": "../../../syntax/StochasticTools/index.html"}>>>]
[]

[Distributions<<<{"href": "../../../syntax/Distributions/index.html"}>>>]
  [R_dist]
    type = Uniform<<<{"description": "Continuous uniform distribution.", "href": "../../../source/distributions/Uniform.html"}>>>
    lower_bound<<<{"description": "Distribution lower bound"}>>> = 1.25E-4
    upper_bound<<<{"description": "Distribution upper bound"}>>> = 1.55E-4
  []
  [power_dist]
    type = Uniform<<<{"description": "Continuous uniform distribution.", "href": "../../../source/distributions/Uniform.html"}>>>
    lower_bound<<<{"description": "Distribution lower bound"}>>> = 60
    upper_bound<<<{"description": "Distribution upper bound"}>>> = 70
  []
[]

[Samplers<<<{"href": "../../../syntax/Samplers/index.html"}>>>]
  [test]
    type = MonteCarlo<<<{"description": "Monte Carlo Sampler.", "href": "../../../source/samplers/MonteCarloSampler.html"}>>>
    num_rows<<<{"description": "The number of rows per matrix to generate."}>>> = 90
    distributions<<<{"description": "The distribution names to be sampled, the number of distributions provided defines the number of columns per matrix."}>>> = 'R_dist power_dist'
    execute_on<<<{"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."}>>> = PRE_MULTIAPP_SETUP
    min_procs_per_row<<<{"description": "This will ensure that the sampler is partitioned properly when 'MultiApp/*/min_procs_per_app' is specified. It is not recommended to use otherwise."}>>> = 2
    max_procs_per_row<<<{"description": "This will ensure that the sampler is partitioned properly when 'MultiApp/*/max_procs_per_app' is specified. It is not recommended to use otherwise."}>>> = 2
    seed<<<{"description": "Random number generator initial seed"}>>>=42
  []
[]

[MultiApps<<<{"href": "../../../syntax/MultiApps/index.html"}>>>]
  [worker]
    type = SamplerFullSolveMultiApp<<<{"description": "Creates a full-solve type sub-application for each row of each Sampler matrix.", "href": "../../../source/multiapps/SamplerFullSolveMultiApp.html"}>>>
    input_files<<<{"description": "The input file for each App.  If this parameter only contains one input file it will be used for all of the Apps.  When using 'positions_from_file' it is also admissable to provide one input_file per file."}>>> = 2d-reconst.i
    sampler<<<{"description": "The Sampler object to utilize for creating MultiApps."}>>> = test
    mode<<<{"description": "The operation mode, 'normal' creates one sub-application for each row in the Sampler and 'batch-reset' and 'batch-restore' creates N sub-applications, where N is the minimum of 'num_rows' in the Sampler and floor(number of processes / min_procs_per_app). To run the rows in the Sampler, 'batch-reset' will destroy and re-create sub-apps as needed, whereas the 'batch-restore' will backup and restore sub-apps to the initial state prior to execution, without destruction."}>>> = batch-reset
    min_procs_per_app<<<{"description": "Minimum number of processors to give to each App in this MultiApp.  Useful for larger, distributed mesh solves."}>>> = 2
    max_procs_per_app<<<{"description": "Maximum number of processors to give to each App in this MultiApp.  Useful for restricting small solves to just a few procs so they don't get spread out"}>>> = 2
  []
[]

[Controls<<<{"href": "../../../syntax/Controls/index.html"}>>>]
  [cmdline]
    type = MultiAppSamplerControl<<<{"description": "Control for modifying the command line arguments of MultiApps.", "href": "../../../source/controls/MultiAppSamplerControl.html"}>>>
    multi_app<<<{"description": "The name of the MultiApp to control."}>>> = worker
    sampler<<<{"description": "The Sampler object to utilize for altering the command line options of the MultiApp."}>>> = test
    param_names<<<{"description": "The names of the command line parameters to set via the sampled data."}>>> = 'R power'
  []
[]

[Transfers<<<{"href": "../../../syntax/Transfers/index.html"}>>>]
  [results]
    type = SamplerReporterTransfer<<<{"description": "Transfers data from Reporters on the sub-application to a StochasticReporter on the main application.", "href": "../../../source/transfers/SamplerReporterTransfer.html"}>>>
    from_multi_app<<<{"description": "The name of the MultiApp to receive data from"}>>> = worker
    sampler<<<{"description": "A the Sampler object that Transfer is associated.."}>>> = test
    stochastic_reporter<<<{"description": "The name of the StochasticReporter object to transfer values to."}>>> = matrix
    from_reporter<<<{"description": "The name(s) of the Reporter(s) on the sub-app to transfer from."}>>> = 'l2error/value'
  []
[]

[Reporters<<<{"href": "../../../syntax/Reporters/index.html"}>>>]
  [matrix]
    type = StochasticMatrix<<<{"description": "Tool for extracting Sampler object data and storing data from stochastic simulations.", "href": "../../../source/reporters/StochasticMatrix.html"}>>>
    sampler<<<{"description": "The sample from which to extract distribution data."}>>> = test
    parallel_type<<<{"description": "This parameter will determine how the stochastic data is gathered. It is common for outputting purposes that this parameter be set to ROOT, otherwise, many files will be produced showing the values on each processor. However, if there are lot of samples, gathering on root may be memory restrictive."}>>> = ROOT
  []
[]

[Outputs<<<{"href": "../../../syntax/Outputs/index.html"}>>>]
  [json]
    type = JSON<<<{"description": "Output for Reporter values using JSON format.", "href": "../../../source/outputs/JSONOutput.html"}>>>
    execute_on<<<{"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."}>>> = FINAL
    execute_system_information_on<<<{"description": "Control when the output of the simulation information occurs"}>>> = NONE
  []
[]

This simply gathers the relative $var element = document.getElementById("moose-equation-7eea7fe2-9b80-4258-a5b9-d820e1d262d6");katex.render("L^2", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ errors from the subapplications and presents the results in a json file. The error histogram is presented below. We see that the maximum integrated $var element = document.getElementById("moose-equation-32575f44-b069-43fa-ae34-e5ffc9d75fef");katex.render("L^2", element, {displayMode:false,throwOnError:false,macros:{"\\eqc":"\\,,","\\eqp":"\\,.","\\pd":"\\frac{\\partial #1}{\\partial #2}","\\pr":"\\left(#1\\right)","\\ddt":"\\frac{d #1}{d t}"}});$ error is below 1.6% with the majority of the samples having error in the 0-0.2% interval.

Figure 5: Histogram of the L2 errors over the 90 test samples.

References

Xuehui Chen, Weihao Mu, Xin Xu, Wei Liu, Lei Huang, and Hao Li. Numerical analysis of double track formation for selective laser melting of 316l stainless steel. Applied Physics A, 127:1–13, 2021.

@article{chen2021numerical,
    author = "Chen, Xuehui and Mu, Weihao and Xu, Xin and Liu, Wei and Huang, Lei and Li, Hao",
    title = "Numerical analysis of double track formation for selective laser melting of 316L stainless steel",
    journal = "Applied Physics A",
    volume = "127",
    pages = "1--13",
    year = "2021",
    publisher = "Springer"
}

Choong S Kim. Thermophysical properties of stainless steels. Technical Report, Argonne National Lab., Ill.(USA), 1975.

@techreport{kim1975thermophysical,
    author = "Kim, Choong S",
    title = "Thermophysical properties of stainless steels",
    year = "1975",
    institution = "Argonne National Lab., Ill.(USA)"
}

David R Noble, Patrick K Notz, Mario J Martinez, and Andrew Michael Kraynik. Use of aria to simulate laser weld pool dynamics for neutron generator production. Technical Report, Sandia National Laboratories (SNL), Albuquerque, NM, and Livermore, CA …, 2007.

@techreport{noble2007use,
    author = "Noble, David R and Notz, Patrick K and Martinez, Mario J and Kraynik, Andrew Michael",
    title = "Use of Aria to simulate laser weld pool dynamics for neutron generator production.",
    year = "2007",
    institution = "Sandia National Laboratories (SNL), Albuquerque, NM, and Livermore, CA\textasciitilde …"
}

Peter Pichler, Thomas Leitner, Erhard Kaschnitz, Johannes Rattenberger, and Gernot Pottlacher. Surface tension and thermal conductivity of nist srm 1155a (aisi 316l stainless steel). International Journal of Thermophysics, 43(5):66, 2022.

@article{pichler2022surface,
    author = "Pichler, Peter and Leitner, Thomas and Kaschnitz, Erhard and Rattenberger, Johannes and Pottlacher, Gernot",
    title = "Surface tension and thermal conductivity of NIST SRM 1155a (AISI 316L stainless steel)",
    journal = "International Journal of Thermophysics",
    volume = "43",
    number = "5",
    pages = "66",
    year = "2022",
    publisher = "Springer"
}

(modules/combined/examples/stochastic/laser_welding_dimred/2d.i)

!include parameters.i

!include mesh.i

[Variables]
  [vel]
    family = LAGRANGE_VEC
  []
  [T]
  []
  [p]
  []
  [disp_x]
  []
  [disp_y]
  []
[]

!include physics_objects.i

[Executioner]
  type = Transient
  end_time = ${endtime}
  dtmin = 1e-10
  dtmax = 1e-5
  petsc_options_iname = '-pc_type -pc_factor_shift_type'
  petsc_options_value = 'lu       NONZERO'
  petsc_options = '-snes_converged_reason -ksp_converged_reason -options_left'
  solve_type = 'NEWTON'
  line_search = 'none'
  nl_max_its = 5
  l_max_its = 100
  [TimeStepper]
    type = IterationAdaptiveDT
    optimal_iterations = 5
    iteration_window = 1
    dt = ${timestep}
    linear_iteration_ratio = 1e6
    growth_factor = 1.25
  []
[]

[Reporters]
  [solution_storage]
    type = SolutionContainer
    execute_on = 'FINAL'
  []
[]

[Debug]
  show_var_residual_norms = true
[]

(modules/combined/examples/stochastic/laser_welding_dimred/parameters.i)

# Process parameters
scanning_speed=1.0 # m/s
power=74 # W (this is the effective power so multiplied by eta)
R=1.25e-4 # m (this is the effective radius)

# Geometric parameters
thickness=1e-4 # m
xmin=-0.1e-3 # m
xmax=0.75e-3 # m
ymin=${fparse -thickness}
surfacetemp=300 # K (temperature at the other side of the plate)

# Time stepping parameters
endtime=4e-4 # s
timestep=${fparse endtime/1000} # s

(modules/combined/examples/stochastic/laser_welding_dimred/mesh.i)

[Mesh]
  [cmg]
    type = GeneratedMeshGenerator
    dim = 2
    xmin = ${xmin}
    xmax = ${xmax}
    ymin = ${fparse ymin}
    ymax = 0
    nx = 161
    ny = 50
  []
  displacements = 'disp_x disp_y'
[]

(modules/combined/examples/stochastic/laser_welding_dimred/physics_objects.i)

[ICs]
  [T]
    type = FunctionIC
    variable = T
    function = '(${surfacetemp} - 300) / ${thickness} * y + ${surfacetemp}'
  []
[]

[Kernels]
  [disp_x]
    type = Diffusion
    variable = disp_x
  []
  [disp_y]
    type = Diffusion
    variable = disp_y
  []
  [mass]
    type = INSADMass
    variable = p
    use_displaced_mesh = true
  []
  [mass_pspg]
    type = INSADMassPSPG
    variable = p
    use_displaced_mesh = true
  []
  [momentum_time]
    type = INSADMomentumTimeDerivative
    variable = vel
    use_displaced_mesh = true
  []
  [momentum_advection]
    type = INSADMomentumAdvection
    variable = vel
    use_displaced_mesh = true
  []
  [momentum_mesh_advection]
    type = INSADMomentumMeshAdvection
    variable = vel
    disp_x = disp_x
    disp_y = disp_y
    use_displaced_mesh = true
  []
  [momentum_viscous]
    type = INSADMomentumViscous
    variable = vel
    use_displaced_mesh = true
  []
  [momentum_pressure]
    type = INSADMomentumPressure
    variable = vel
    pressure = p
    integrate_p_by_parts = true
    use_displaced_mesh = true
  []
  [momentum_supg]
    type = INSADMomentumSUPG
    variable = vel
    material_velocity = relative_velocity
    use_displaced_mesh = true
  []
  [temperature_time]
    type = INSADHeatConductionTimeDerivative
    variable = T
    use_displaced_mesh = true
  []
  [temperature_advection]
    type = INSADEnergyAdvection
    variable = T
    use_displaced_mesh = true
  []
  [temperature_mesh_advection]
    type = INSADEnergyMeshAdvection
    variable = T
    disp_x = disp_x
    disp_y = disp_y
    use_displaced_mesh = true
  []
  [temperature_conduction]
    type = ADHeatConduction
    variable = T
    thermal_conductivity = 'k'
    use_displaced_mesh = true
  []
  [temperature_supg]
    type = INSADEnergySUPG
    variable = T
    velocity = vel
    use_displaced_mesh = true
  []
[]

[BCs]
  [x_no_disp]
    type = DirichletBC
    variable = disp_x
    boundary = 'bottom'
    value = 0
  []
  [y_no_disp]
    type = DirichletBC
    variable = disp_y
    boundary = 'bottom'
    value = 0
  []
  [no_slip]
    type = ADVectorFunctionDirichletBC
    variable = vel
    boundary = 'bottom right left'
  []
  [T_cold]
    type = DirichletBC
    variable = T
    boundary = 'bottom'
    value = 300
  []
  [radiation_flux]
    type = FunctionRadiativeBC
    variable = T
    boundary = 'top'
    emissivity_function = '1'
    Tinfinity = 300
    stefan_boltzmann_constant = 5.67e-8
    use_displaced_mesh = true
  []
  [weld_flux]
    type = GaussianEnergyFluxBC
    variable = T
    boundary = 'top'
    P0 = ${power}
    R = ${R}
    x_beam_coord = '${scanning_speed}*t'
    y_beam_coord = '0'
    use_displaced_mesh = true
  []
  [vapor_recoil]
    type = INSADVaporRecoilPressureMomentumFluxBC
    variable = vel
    boundary = 'top'
    use_displaced_mesh = true
  []
  [surface_tension]
    type = INSADSurfaceTensionBC
    variable = vel
    boundary = 'top'
    use_displaced_mesh = true
  []
  [displace_x_top]
    type = INSADDisplaceBoundaryBC
    boundary = 'top'
    variable = 'disp_x'
    velocity = 'vel'
    component = 0
    associated_subdomain = 0
  []
  [displace_y_top]
    type = INSADDisplaceBoundaryBC
    boundary = 'top'
    variable = 'disp_y'
    velocity = 'vel'
    component = 1
    associated_subdomain = 0
  []
  [displace_x_top_dummy]
    type = INSADDummyDisplaceBoundaryIntegratedBC
    boundary = 'top'
    variable = 'disp_x'
    velocity = 'vel'
    component = 0
  []
  [displace_y_top_dummy]
    type = INSADDummyDisplaceBoundaryIntegratedBC
    boundary = 'top'
    variable = 'disp_y'
    velocity = 'vel'
    component = 1
  []
[]

[Materials]
  [ins_mat]
    type = INSADStabilized3Eqn
    velocity = vel
    pressure = p
    temperature = T
    use_displaced_mesh = true
  []
  [steel]
    type = LaserWeld316LStainlessSteel
    temperature = T
    use_displaced_mesh = true
  []
  [steel_boundary]
    type = LaserWeld316LStainlessSteelBoundary
    boundary = 'top'
    temperature = T
    use_displaced_mesh = true
  []
  [const]
    type = GenericConstantMaterial
    prop_names = 'abs sb_constant'
    prop_values = '1 5.67e-8'
    use_displaced_mesh = true
  []
[]

(modules/combined/examples/stochastic/laser_welding_dimred/train.i)

[StochasticTools]
[]

[Distributions]
  [R_dist]
    type = Uniform
    lower_bound = 1.25E-4
    upper_bound = 1.55E-4
  []
  [power_dist]
    type = Uniform
    lower_bound = 60
    upper_bound = 74
  []
[]

[Samplers]
  [train]
    type = MonteCarlo
    num_rows = 45
    distributions = 'R_dist power_dist'
    execute_on = PRE_MULTIAPP_SETUP
    min_procs_per_row = 2
    max_procs_per_row = 2
  []
[]

[MultiApps]
  [worker]
    type = SamplerFullSolveMultiApp
    input_files = 2d.i
    sampler = train
    mode = batch-reset
    min_procs_per_app = 2
    max_procs_per_app = 2
  []
[]

[Controls]
  [cmdline]
    type = MultiAppSamplerControl
    multi_app = worker
    sampler = train
    param_names = 'R power'
  []
[]

[Transfers]
  [solution_transfer]
    type = SerializedSolutionTransfer
    parallel_storage = parallel_storage
    from_multi_app = worker
    sampler = train
    solution_container = solution_storage
    variables = "T"
    serialize_on_root = true
  []
[]

[Reporters]
  [parallel_storage]
    type = ParallelSolutionStorage
    variables = "T"
    outputs = none
  []
  [reduced_sol]
    type = MappingReporter
    sampler = train
    parallel_storage = parallel_storage
    mapping = pod
    variables = "T"
    outputs = json
    execute_on = final
    parallel_type = ROOT
  []
  [matrix]
    type = StochasticMatrix
    sampler = train
    parallel_type = ROOT
  []
  [svd]
    type = SingularTripletReporter
    pod_mapping = pod
    variables = "T"
    execute_on = final
  []
[]

[VariableMappings]
  [pod]
    type = PODMapping
    solution_storage = parallel_storage
    variables = "T"
    num_modes_to_compute = '8'
    extra_slepc_options = "-svd_monitor_all"
  []
[]

[Trainers]
  [mogp]
    type = GaussianProcessTrainer
    execute_on = final
    covariance_function = 'lmc'
    standardize_params = 'true'
    standardize_data = 'true'
    sampler = train
    response_type = vector_real
    response = reduced_sol/T_pod
    tune_parameters = 'lmc:acoeff_0 lmc:lambda_0 covar:signal_variance covar:length_factor'
    tuning_min = '1e-9 1e-9 1e-9 1e-9'
    tuning_max = '1e16 1e16 1e16  1e16'
    num_iters = 10000
    learning_rate = 0.0005
    show_every_nth_iteration = 10
  []
[]

[Covariance]
  [covar]
    type = SquaredExponentialCovariance
    signal_variance = 1.0
    noise_variance = 0.0
    length_factor = '0.1 0.1'
  []
  [lmc]
    type = LMC
    covariance_functions = covar
    num_outputs = 8
    num_latent_funcs = 1
  []
[]


[Outputs]
  [json]
    type = JSON
    execute_on = FINAL
    execute_system_information_on = NONE
  []
  [pod_out]
    type = MappingOutput
    mappings = pod
    execute_on = FINAL
  []
  [mogp_out]
    type = SurrogateTrainerOutput
    trainers = mogp
    execute_on = FINAL
  []
[]

(modules/combined/examples/stochastic/laser_welding_dimred/train.i)

[StochasticTools]
[]

[Distributions]
  [R_dist]
    type = Uniform
    lower_bound = 1.25E-4
    upper_bound = 1.55E-4
  []
  [power_dist]
    type = Uniform
    lower_bound = 60
    upper_bound = 74
  []
[]

[Samplers]
  [train]
    type = MonteCarlo
    num_rows = 45
    distributions = 'R_dist power_dist'
    execute_on = PRE_MULTIAPP_SETUP
    min_procs_per_row = 2
    max_procs_per_row = 2
  []
[]

[MultiApps]
  [worker]
    type = SamplerFullSolveMultiApp
    input_files = 2d.i
    sampler = train
    mode = batch-reset
    min_procs_per_app = 2
    max_procs_per_app = 2
  []
[]

[Controls]
  [cmdline]
    type = MultiAppSamplerControl
    multi_app = worker
    sampler = train
    param_names = 'R power'
  []
[]

[Transfers]
  [solution_transfer]
    type = SerializedSolutionTransfer
    parallel_storage = parallel_storage
    from_multi_app = worker
    sampler = train
    solution_container = solution_storage
    variables = "T"
    serialize_on_root = true
  []
[]

[Reporters]
  [parallel_storage]
    type = ParallelSolutionStorage
    variables = "T"
    outputs = none
  []
  [reduced_sol]
    type = MappingReporter
    sampler = train
    parallel_storage = parallel_storage
    mapping = pod
    variables = "T"
    outputs = json
    execute_on = final
    parallel_type = ROOT
  []
  [matrix]
    type = StochasticMatrix
    sampler = train
    parallel_type = ROOT
  []
  [svd]
    type = SingularTripletReporter
    pod_mapping = pod
    variables = "T"
    execute_on = final
  []
[]

[VariableMappings]
  [pod]
    type = PODMapping
    solution_storage = parallel_storage
    variables = "T"
    num_modes_to_compute = '8'
    extra_slepc_options = "-svd_monitor_all"
  []
[]

[Trainers]
  [mogp]
    type = GaussianProcessTrainer
    execute_on = final
    covariance_function = 'lmc'
    standardize_params = 'true'
    standardize_data = 'true'
    sampler = train
    response_type = vector_real
    response = reduced_sol/T_pod
    tune_parameters = 'lmc:acoeff_0 lmc:lambda_0 covar:signal_variance covar:length_factor'
    tuning_min = '1e-9 1e-9 1e-9 1e-9'
    tuning_max = '1e16 1e16 1e16  1e16'
    num_iters = 10000
    learning_rate = 0.0005
    show_every_nth_iteration = 10
  []
[]

[Covariance]
  [covar]
    type = SquaredExponentialCovariance
    signal_variance = 1.0
    noise_variance = 0.0
    length_factor = '0.1 0.1'
  []
  [lmc]
    type = LMC
    covariance_functions = covar
    num_outputs = 8
    num_latent_funcs = 1
  []
[]


[Outputs]
  [json]
    type = JSON
    execute_on = FINAL
    execute_system_information_on = NONE
  []
  [pod_out]
    type = MappingOutput
    mappings = pod
    execute_on = FINAL
  []
  [mogp_out]
    type = SurrogateTrainerOutput
    trainers = mogp
    execute_on = FINAL
  []
[]

(modules/combined/examples/stochastic/laser_welding_dimred/train.i)

[StochasticTools]
[]

[Distributions]
  [R_dist]
    type = Uniform
    lower_bound = 1.25E-4
    upper_bound = 1.55E-4
  []
  [power_dist]
    type = Uniform
    lower_bound = 60
    upper_bound = 74
  []
[]

[Samplers]
  [train]
    type = MonteCarlo
    num_rows = 45
    distributions = 'R_dist power_dist'
    execute_on = PRE_MULTIAPP_SETUP
    min_procs_per_row = 2
    max_procs_per_row = 2
  []
[]

[MultiApps]
  [worker]
    type = SamplerFullSolveMultiApp
    input_files = 2d.i
    sampler = train
    mode = batch-reset
    min_procs_per_app = 2
    max_procs_per_app = 2
  []
[]

[Controls]
  [cmdline]
    type = MultiAppSamplerControl
    multi_app = worker
    sampler = train
    param_names = 'R power'
  []
[]

[Transfers]
  [solution_transfer]
    type = SerializedSolutionTransfer
    parallel_storage = parallel_storage
    from_multi_app = worker
    sampler = train
    solution_container = solution_storage
    variables = "T"
    serialize_on_root = true
  []
[]

[Reporters]
  [parallel_storage]
    type = ParallelSolutionStorage
    variables = "T"
    outputs = none
  []
  [reduced_sol]
    type = MappingReporter
    sampler = train
    parallel_storage = parallel_storage
    mapping = pod
    variables = "T"
    outputs = json
    execute_on = final
    parallel_type = ROOT
  []
  [matrix]
    type = StochasticMatrix
    sampler = train
    parallel_type = ROOT
  []
  [svd]
    type = SingularTripletReporter
    pod_mapping = pod
    variables = "T"
    execute_on = final
  []
[]

[VariableMappings]
  [pod]
    type = PODMapping
    solution_storage = parallel_storage
    variables = "T"
    num_modes_to_compute = '8'
    extra_slepc_options = "-svd_monitor_all"
  []
[]

[Trainers]
  [mogp]
    type = GaussianProcessTrainer
    execute_on = final
    covariance_function = 'lmc'
    standardize_params = 'true'
    standardize_data = 'true'
    sampler = train
    response_type = vector_real
    response = reduced_sol/T_pod
    tune_parameters = 'lmc:acoeff_0 lmc:lambda_0 covar:signal_variance covar:length_factor'
    tuning_min = '1e-9 1e-9 1e-9 1e-9'
    tuning_max = '1e16 1e16 1e16  1e16'
    num_iters = 10000
    learning_rate = 0.0005
    show_every_nth_iteration = 10
  []
[]

[Covariance]
  [covar]
    type = SquaredExponentialCovariance
    signal_variance = 1.0
    noise_variance = 0.0
    length_factor = '0.1 0.1'
  []
  [lmc]
    type = LMC
    covariance_functions = covar
    num_outputs = 8
    num_latent_funcs = 1
  []
[]


[Outputs]
  [json]
    type = JSON
    execute_on = FINAL
    execute_system_information_on = NONE
  []
  [pod_out]
    type = MappingOutput
    mappings = pod
    execute_on = FINAL
  []
  [mogp_out]
    type = SurrogateTrainerOutput
    trainers = mogp
    execute_on = FINAL
  []
[]

(modules/combined/examples/stochastic/laser_welding_dimred/train.i)

[StochasticTools]
[]

[Distributions]
  [R_dist]
    type = Uniform
    lower_bound = 1.25E-4
    upper_bound = 1.55E-4
  []
  [power_dist]
    type = Uniform
    lower_bound = 60
    upper_bound = 74
  []
[]

[Samplers]
  [train]
    type = MonteCarlo
    num_rows = 45
    distributions = 'R_dist power_dist'
    execute_on = PRE_MULTIAPP_SETUP
    min_procs_per_row = 2
    max_procs_per_row = 2
  []
[]

[MultiApps]
  [worker]
    type = SamplerFullSolveMultiApp
    input_files = 2d.i
    sampler = train
    mode = batch-reset
    min_procs_per_app = 2
    max_procs_per_app = 2
  []
[]

[Controls]
  [cmdline]
    type = MultiAppSamplerControl
    multi_app = worker
    sampler = train
    param_names = 'R power'
  []
[]

[Transfers]
  [solution_transfer]
    type = SerializedSolutionTransfer
    parallel_storage = parallel_storage
    from_multi_app = worker
    sampler = train
    solution_container = solution_storage
    variables = "T"
    serialize_on_root = true
  []
[]

[Reporters]
  [parallel_storage]
    type = ParallelSolutionStorage
    variables = "T"
    outputs = none
  []
  [reduced_sol]
    type = MappingReporter
    sampler = train
    parallel_storage = parallel_storage
    mapping = pod
    variables = "T"
    outputs = json
    execute_on = final
    parallel_type = ROOT
  []
  [matrix]
    type = StochasticMatrix
    sampler = train
    parallel_type = ROOT
  []
  [svd]
    type = SingularTripletReporter
    pod_mapping = pod
    variables = "T"
    execute_on = final
  []
[]

[VariableMappings]
  [pod]
    type = PODMapping
    solution_storage = parallel_storage
    variables = "T"
    num_modes_to_compute = '8'
    extra_slepc_options = "-svd_monitor_all"
  []
[]

[Trainers]
  [mogp]
    type = GaussianProcessTrainer
    execute_on = final
    covariance_function = 'lmc'
    standardize_params = 'true'
    standardize_data = 'true'
    sampler = train
    response_type = vector_real
    response = reduced_sol/T_pod
    tune_parameters = 'lmc:acoeff_0 lmc:lambda_0 covar:signal_variance covar:length_factor'
    tuning_min = '1e-9 1e-9 1e-9 1e-9'
    tuning_max = '1e16 1e16 1e16  1e16'
    num_iters = 10000
    learning_rate = 0.0005
    show_every_nth_iteration = 10
  []
[]

[Covariance]
  [covar]
    type = SquaredExponentialCovariance
    signal_variance = 1.0
    noise_variance = 0.0
    length_factor = '0.1 0.1'
  []
  [lmc]
    type = LMC
    covariance_functions = covar
    num_outputs = 8
    num_latent_funcs = 1
  []
[]


[Outputs]
  [json]
    type = JSON
    execute_on = FINAL
    execute_system_information_on = NONE
  []
  [pod_out]
    type = MappingOutput
    mappings = pod
    execute_on = FINAL
  []
  [mogp_out]
    type = SurrogateTrainerOutput
    trainers = mogp
    execute_on = FINAL
  []
[]

(modules/combined/examples/stochastic/laser_welding_dimred/2d-reconst.i)

!include parameters.i

!include mesh.i

[Variables]
  [vel]
    family = LAGRANGE_VEC
  []
  [T]
  []
  [p]
  []
  [disp_x]
  []
  [disp_y]
  []
[]

[AuxVariables]
  [T_reconst]
  []
[]

!include physics_objects.i

[UserObjects]
  [inverse]
    type = InverseMapping
    mapping = pod
    variable_to_fill = "T_reconst"
    variable_to_reconstruct = "T"
    surrogate = mogp
    parameters = '${R} ${power}'
    execute_on = TIMESTEP_END
  []
[]

[Surrogates]
  [mogp]
    type = GaussianProcessSurrogate
    filename = "train_mogp_out_mogp.rd"
  []
[]

[VariableMappings]
  [pod]
    type = PODMapping
    filename = "train_pod_out_pod.rd"
  []
[]

[Executioner]
  type = Transient
  end_time = ${endtime}
  dtmin = 1e-10
  dtmax = 1e-5
  petsc_options_iname = '-pc_type -pc_factor_shift_type'
  petsc_options_value = 'lu       NONZERO'
  petsc_options = '-snes_converged_reason -ksp_converged_reason -options_left'
  solve_type = 'NEWTON'
  line_search = 'none'
  nl_max_its = 16
  l_max_its = 100
  [TimeStepper]
    type = IterationAdaptiveDT
    optimal_iterations = 10
    iteration_window = 4
    dt = ${timestep}
    linear_iteration_ratio = 1e6
    growth_factor = 1.25
  []
[]

[Postprocessors]
  [l2error]
    type = ElementL2Difference
    variable = T
    other_variable = T_reconst
  []
[]

[Debug]
  show_var_residual_norms = true
[]

[Postprocessors]
[]

(modules/combined/examples/stochastic/laser_welding_dimred/2d-reconst.i)

!include parameters.i

!include mesh.i

[Variables]
  [vel]
    family = LAGRANGE_VEC
  []
  [T]
  []
  [p]
  []
  [disp_x]
  []
  [disp_y]
  []
[]

[AuxVariables]
  [T_reconst]
  []
[]

!include physics_objects.i

[UserObjects]
  [inverse]
    type = InverseMapping
    mapping = pod
    variable_to_fill = "T_reconst"
    variable_to_reconstruct = "T"
    surrogate = mogp
    parameters = '${R} ${power}'
    execute_on = TIMESTEP_END
  []
[]

[Surrogates]
  [mogp]
    type = GaussianProcessSurrogate
    filename = "train_mogp_out_mogp.rd"
  []
[]

[VariableMappings]
  [pod]
    type = PODMapping
    filename = "train_pod_out_pod.rd"
  []
[]

[Executioner]
  type = Transient
  end_time = ${endtime}
  dtmin = 1e-10
  dtmax = 1e-5
  petsc_options_iname = '-pc_type -pc_factor_shift_type'
  petsc_options_value = 'lu       NONZERO'
  petsc_options = '-snes_converged_reason -ksp_converged_reason -options_left'
  solve_type = 'NEWTON'
  line_search = 'none'
  nl_max_its = 16
  l_max_its = 100
  [TimeStepper]
    type = IterationAdaptiveDT
    optimal_iterations = 10
    iteration_window = 4
    dt = ${timestep}
    linear_iteration_ratio = 1e6
    growth_factor = 1.25
  []
[]

[Postprocessors]
  [l2error]
    type = ElementL2Difference
    variable = T
    other_variable = T_reconst
  []
[]

[Debug]
  show_var_residual_norms = true
[]

[Postprocessors]
[]

(modules/combined/examples/stochastic/laser_welding_dimred/2d-reconst.i)

!include parameters.i

!include mesh.i

[Variables]
  [vel]
    family = LAGRANGE_VEC
  []
  [T]
  []
  [p]
  []
  [disp_x]
  []
  [disp_y]
  []
[]

[AuxVariables]
  [T_reconst]
  []
[]

!include physics_objects.i

[UserObjects]
  [inverse]
    type = InverseMapping
    mapping = pod
    variable_to_fill = "T_reconst"
    variable_to_reconstruct = "T"
    surrogate = mogp
    parameters = '${R} ${power}'
    execute_on = TIMESTEP_END
  []
[]

[Surrogates]
  [mogp]
    type = GaussianProcessSurrogate
    filename = "train_mogp_out_mogp.rd"
  []
[]

[VariableMappings]
  [pod]
    type = PODMapping
    filename = "train_pod_out_pod.rd"
  []
[]

[Executioner]
  type = Transient
  end_time = ${endtime}
  dtmin = 1e-10
  dtmax = 1e-5
  petsc_options_iname = '-pc_type -pc_factor_shift_type'
  petsc_options_value = 'lu       NONZERO'
  petsc_options = '-snes_converged_reason -ksp_converged_reason -options_left'
  solve_type = 'NEWTON'
  line_search = 'none'
  nl_max_its = 16
  l_max_its = 100
  [TimeStepper]
    type = IterationAdaptiveDT
    optimal_iterations = 10
    iteration_window = 4
    dt = ${timestep}
    linear_iteration_ratio = 1e6
    growth_factor = 1.25
  []
[]

[Postprocessors]
  [l2error]
    type = ElementL2Difference
    variable = T
    other_variable = T_reconst
  []
[]

[Debug]
  show_var_residual_norms = true
[]

[Postprocessors]
[]

(modules/combined/examples/stochastic/laser_welding_dimred/2d-reconst.i)

!include parameters.i

!include mesh.i

[Variables]
  [vel]
    family = LAGRANGE_VEC
  []
  [T]
  []
  [p]
  []
  [disp_x]
  []
  [disp_y]
  []
[]

[AuxVariables]
  [T_reconst]
  []
[]

!include physics_objects.i

[UserObjects]
  [inverse]
    type = InverseMapping
    mapping = pod
    variable_to_fill = "T_reconst"
    variable_to_reconstruct = "T"
    surrogate = mogp
    parameters = '${R} ${power}'
    execute_on = TIMESTEP_END
  []
[]

[Surrogates]
  [mogp]
    type = GaussianProcessSurrogate
    filename = "train_mogp_out_mogp.rd"
  []
[]

[VariableMappings]
  [pod]
    type = PODMapping
    filename = "train_pod_out_pod.rd"
  []
[]

[Executioner]
  type = Transient
  end_time = ${endtime}
  dtmin = 1e-10
  dtmax = 1e-5
  petsc_options_iname = '-pc_type -pc_factor_shift_type'
  petsc_options_value = 'lu       NONZERO'
  petsc_options = '-snes_converged_reason -ksp_converged_reason -options_left'
  solve_type = 'NEWTON'
  line_search = 'none'
  nl_max_its = 16
  l_max_its = 100
  [TimeStepper]
    type = IterationAdaptiveDT
    optimal_iterations = 10
    iteration_window = 4
    dt = ${timestep}
    linear_iteration_ratio = 1e6
    growth_factor = 1.25
  []
[]

[Postprocessors]
  [l2error]
    type = ElementL2Difference
    variable = T
    other_variable = T_reconst
  []
[]

[Debug]
  show_var_residual_norms = true
[]

[Postprocessors]
[]

(modules/combined/examples/stochastic/laser_welding_dimred/test.i)

[StochasticTools]
[]

[Distributions]
  [R_dist]
    type = Uniform
    lower_bound = 1.25E-4
    upper_bound = 1.55E-4
  []
  [power_dist]
    type = Uniform
    lower_bound = 60
    upper_bound = 70
  []
[]

[Samplers]
  [test]
    type = MonteCarlo
    num_rows = 90
    distributions = 'R_dist power_dist'
    execute_on = PRE_MULTIAPP_SETUP
    min_procs_per_row = 2
    max_procs_per_row = 2
    seed=42
  []
[]

[MultiApps]
  [worker]
    type = SamplerFullSolveMultiApp
    input_files = 2d-reconst.i
    sampler = test
    mode = batch-reset
    min_procs_per_app = 2
    max_procs_per_app = 2
  []
[]

[Controls]
  [cmdline]
    type = MultiAppSamplerControl
    multi_app = worker
    sampler = test
    param_names = 'R power'
  []
[]

[Transfers]
  [results]
    type = SamplerReporterTransfer
    from_multi_app = worker
    sampler = test
    stochastic_reporter = matrix
    from_reporter = 'l2error/value'
  []
[]

[Reporters]
  [matrix]
    type = StochasticMatrix
    sampler = test
    parallel_type = ROOT
  []
[]

[Outputs]
  [json]
    type = JSON
    execute_on = FINAL
    execute_system_information_on = NONE
  []
[]

Problem Description
Results
Reference results
Training the reduced - order model
Evaluating the error of the ROM
References