Force Inversion Example: Neumann Boundary Condition

Background

The MOOSE optimization module provides a flexible framework for solving inverse optimization problems in MOOSE. This page is part of a set of examples for different types of inverse optimization problems.

Example: Neumann Boundary Condition

Main-App Optimization Executioner

In this example, we will be optimizing the multi-objective function given by:

(1)

where represents the displacement vector, denotes the parameters, and are the observed and target -components of the displacement at the -th point, and and are the observed and target -components of the displacement at the -th point.

The derivative of the objective function with respect to the parameters can be expressed as:

(2)

In this example we will be optimizing the two forcing parameters on the right side of a cantilever beam. The forcing parameters control the force that is applied to the beam in each direction.

The main app, shown in Listing 1, only receives the final objective and gradient.

Listing 1: MainApp Input File

[Optimization]
[]

[Mesh]
  type = GeneratedMesh
  dim = 2
  nx = 10
  ny = 2
  xmin = 0.0
  xmax = 5.0
  ymin = 0.0
  ymax = 1.0
[]

[OptimizationReporter]
  type = GeneralOptimization
  objective_name = objective_value
  parameter_names = 'force_right'
  num_values = '2'
  initial_condition = '100'
  outputs = 'csv'
[]

[Executioner]
  type = Optimize
  tao_solver = taobqnktr
  petsc_options_iname = '-tao_gttol -tao_max_it'
  petsc_options_value = '1e-7 500 '
  verbose = true
[]

[MultiApps]
  [forward]
    type = FullSolveMultiApp
    input_files = forward_and_adjoint.i
    execute_on = "FORWARD"
    clone_parent_mesh = true
  []
[]

[Transfers]
  [toForward_measument]
    type = MultiAppReporterTransfer
    to_multi_app = forward
    from_reporters = 'OptimizationReporter/force_right'
    to_reporters = 'params/right_values'
  []
  [fromForward]
    type = MultiAppReporterTransfer
    from_multi_app = forward
    from_reporters = 'combined/gradient obj/obj_val'
    to_reporters = 'OptimizationReporter/grad_force_right OptimizationReporter/objective_value'
  []
[]

[Outputs]
  csv = true
  execute_on = 'FINAL'
[]
(modules/combined/test/tests/optimization/invOpt_mechanics/main_auto_adjoint.i)

Forward and Adjoint Problem Sub-App

The definition of the objective function in Eq. (1) allows for the gradient to be the addition of two gradients: one for the -component and one for the -component. The forward sub app is responsible for correctly creating both sub-gradients and combining them, as well as the objective functions, to be used by the main app. This example will only go into detail on how to calculate and combine these objectives and gradients. For more details on setting up the entire optimization process, look at Example 1: Point Loads.

For functions, we split the parameter vector, which has two entries, into two separate functions: one that controls the force for the -component and another that controls the force for the -component.

Listing 2: Functions

[Functions]
  [right_fy_func]
    type = ParsedOptimizationFunction
    expression = 'val_y'
    param_symbol_names = 'val_x val_y'
    param_vector_name = 'params/right_values'
  []
  [right_fx_func]
    type = ParsedOptimizationFunction
    expression = 'val_x'
    param_symbol_names = 'val_x val_y'
    param_vector_name = 'params/right_values'
  []
[]
(modules/combined/test/tests/optimization/invOpt_mechanics/forward_and_adjoint.i)

For the boundary conditions, we apply the forces as Neumann boundary conditions.

Listing 3: Boundary Conditions

[BCs]
  [left_ux]
    type = DirichletBC
    variable = disp_x
    boundary = left
    value = 0
  []
  [left_uy]
    type = DirichletBC
    variable = disp_y
    boundary = left
    value = 0
  []
  [right_fy]
    type = FunctionNeumannBC
    variable = disp_y
    boundary = right
    function = right_fy_func
  []
  [right_fx]
    type = FunctionNeumannBC
    variable = disp_x
    boundary = right
    function = right_fx_func
  []
[]
(modules/combined/test/tests/optimization/invOpt_mechanics/forward_and_adjoint.i)

For the reporters, we create two OptimizationData objects that hold the observed measurement data for and individually. This Reporter will calculate the individual objective values. The locations or the number of points for and don't need to match. Also, we use a ParsedVectorReporter and a ParsedScalarReporter to combine the individual gradients and objective function, respectively. These are the final reporters that will be transferred to the main app.

Listing 4: Reporters

[Reporters]
  [measure_data_x]
    type = OptimizationData
    objective_name = objective_value
    variable = disp_x
    measurement_points = '5.0 1.0 0.0'
    measurement_values = '-13.00873469088'
  []
  [measure_data_y]
    type = OptimizationData
    objective_name = objective_value
    variable = disp_y
    measurement_points = '5.0 1.0 0.0'
    measurement_values = '85.008027719915'
  []
  [params]
    type = ConstantReporter
    real_vector_names = 'right_values'
    real_vector_values = '-1300 2100 ' # True Values
  []
  [combined]
    type = ParsedVectorReporter
    name = gradient
    reporter_names = 'adjoint_pt_x/inner_product adjoint_pt_y/inner_product'
    reporter_symbols = 'a b'
    expression = 'a+b'
    execute_on = ADJOINT_TIMESTEP_END
    # Just to confirm this happens after the gradient calcutions
    execution_order_group = 1
  []
  [obj]
    type = ParsedScalarReporter
    name = obj_val
    reporter_names = 'measure_data_x/objective_value
    measure_data_y/objective_value'
    reporter_symbols = 'a b'
    expression = 'a+b'
    execute_on = ADJOINT_TIMESTEP_END
    # Just to confirm this happens after the gradient calcutions
    execution_order_group = 1
  []
[]
(modules/combined/test/tests/optimization/invOpt_mechanics/forward_and_adjoint.i)

The Dirac kernels show how we apply the misfit to the adjoint variables, with the adjoint getting the misfit from the data and the getting the misfit from the data.

Listing 5: Dirac Kernels

[DiracKernels]
  [pt_x]
    type = ReporterPointSource
    variable = adjoint_x
    x_coord_name = measure_data_x/measurement_xcoord
    y_coord_name = measure_data_x/measurement_ycoord
    z_coord_name = measure_data_x/measurement_zcoord
    value_name = measure_data_x/misfit_values
  []
  [pt_y]
    type = ReporterPointSource
    variable = adjoint_y
    x_coord_name = measure_data_y/measurement_xcoord
    y_coord_name = measure_data_y/measurement_ycoord
    z_coord_name = measure_data_y/measurement_zcoord
    value_name = measure_data_y/misfit_values
  []
[]
(modules/combined/test/tests/optimization/invOpt_mechanics/forward_and_adjoint.i)

Finally, in vector postprocessors, we show that we calculate the gradients for both components using SideOptimizationNeumannFunctionInnerProduct with the correct adjoint variable and function. These gradients are then combined in the reporters block as previously discussed.

Listing 6: Vector Postprocessors

[VectorPostprocessors]
  [point_sample]
    type = PointValueSampler
    variable = 'disp_x disp_y'
    points = '5.0 1.0 0'
    sort_by = x
  []

  [adjoint_pt_x]
    type = SideOptimizationNeumannFunctionInnerProduct
    variable = adjoint_x
    function = right_fx_func
    boundary = right
    execute_on = ADJOINT_TIMESTEP_END
  []
  [adjoint_pt_y]
    type = SideOptimizationNeumannFunctionInnerProduct
    variable = adjoint_y
    function = right_fy_func
    boundary = right
    execute_on = ADJOINT_TIMESTEP_END
  []
[]
(modules/combined/test/tests/optimization/invOpt_mechanics/forward_and_adjoint.i)