Skip to content

Commit

Permalink
bimaterial inversion using adjoint based gradient closes idaholab#30
Browse files Browse the repository at this point in the history
  • Loading branch information
aaelmeli authored and dschwen committed Sep 13, 2022
1 parent 0bfd709 commit 988ae03
Show file tree
Hide file tree
Showing 5 changed files with 591 additions and 0 deletions.
Original file line number Diff line number Diff line change
@@ -0,0 +1,243 @@
[Mesh]
[gmg]
type = GeneratedMeshGenerator
dim = 2
nx = 16
ny = 16
xmin = -4
xmax = 4
ymin = -4
ymax = 4
[]
[bimaterial]
type = SubdomainBoundingBoxGenerator
input = gmg
block_id = 1
bottom_left = '-100 -100 -100'
top_right = '100 0 100'
[]
[name_blocks]
type = RenameBlockGenerator
input = bimaterial
old_block = '0 1'
new_block = 'top bottom'
[]
[]

[Variables]
[temperature]
[]
[]

[Kernels]
[conduction]
type = MatDiffusion
diffusivity = diffusivity
variable = temperature
[]
[]


[Reporters]
[misfit]
type=OptimizationData
[]
[]

[DiracKernels]
[pt]
type = ReporterPointSource
variable = temperature
x_coord_name = misfit/measurement_xcoord
y_coord_name = misfit/measurement_ycoord
z_coord_name = misfit/measurement_zcoord
value_name = misfit/misfit_values
[]
[]

[BCs]
[bottom]
type = DirichletBC
variable = temperature
boundary = bottom
value = 0
[]
[]

[AuxVariables]
[forwardAdjoint]
[]
[temperature_forward]
[]
[grad_Tx]
order = CONSTANT
family = MONOMIAL
[]
[grad_Ty]
order = CONSTANT
family = MONOMIAL
[]
[grad_Tz]
order = CONSTANT
family = MONOMIAL
[]
[grad_Tfx]
order = CONSTANT
family = MONOMIAL
[]
[grad_Tfy]
order = CONSTANT
family = MONOMIAL
[]
[grad_Tfz]
order = CONSTANT
family = MONOMIAL
[]
[gradient]
order = CONSTANT
family = MONOMIAL
[]
[]

[AuxKernels]
[grad_Tx]
type = VariableGradientComponent
component = x
variable = grad_Tx
gradient_variable = temperature
[]
[grad_Ty]
type = VariableGradientComponent
component = y
variable = grad_Ty
gradient_variable = temperature
[]
[grad_Tz]
type = VariableGradientComponent
component = z
variable = grad_Tz
gradient_variable = temperature
[]
[grad_Tfx]
type = VariableGradientComponent
component = x
variable = grad_Tfx
gradient_variable = temperature_forward
[]
[grad_Tfy]
type = VariableGradientComponent
component = y
variable = grad_Tfy
gradient_variable = temperature_forward
[]
[grad_Tfz]
type = VariableGradientComponent
component = z
variable = grad_Tfz
gradient_variable = temperature_forward
[]
[gradient]
type = ParsedAux
variable = gradient
args = 'grad_Tx grad_Ty grad_Tz grad_Tfx grad_Tfy grad_Tfz'
function = '-grad_Tx*grad_Tfx-grad_Ty*grad_Tfy-grad_Tz*grad_Tfz'#we need to include the material derivative, which can be captured when computing the flux based on the derivative of the material.
[]
[forwardAdjoint] # I am not sure why do we need this?
type = ParsedAux
variable = forwardAdjoint
args = 'temperature_forward temperature'
function = 'temperature_forward*temperature'
[]
[]

[Functions]
[diffusivity_top_function]
type = ParsedFunction
value = alpha
vars = alpha
vals = d_top
[]
[diffusivity_bottom_function]
type = ParsedFunction
value = alpha
vars = alpha
vals = d_bot
[]
[]

[Materials] #same material as what was used in the forward model
[mat_top]
type = GenericFunctionMaterial
block = top
prop_names = diffusivity
prop_values = diffusivity_top_function
[]
[mat_bottom]
type = GenericFunctionMaterial
block = bottom
prop_names = diffusivity
prop_values = diffusivity_bottom_function
[]
[]

[Postprocessors]
[d_bot]
type = VectorPostprocessorComponent
index = 0
vectorpostprocessor = vector_pp
vector_name = diffusivity_values
execute_on = 'linear'
[]
[d_top]
type = VectorPostprocessorComponent
index = 1
vectorpostprocessor = vector_pp
vector_name = diffusivity_values
execute_on = 'linear'
[]
############
# we need to combine the two in one vector.
[grad_bottom] #compute the integral of the gradient variable on the bottom block (first parameter)
type = ElementIntegralVariablePostprocessor
variable = gradient
execute_on = 'final'
block=bottom
outputs=grad_bottom
[]
[grad_top] #compute the integral of the gradient variable on the bottom block (second parameter)
type = ElementIntegralVariablePostprocessor
variable = gradient
execute_on = 'final'
block=top
outputs=grad_top
[]
############
[]

[VectorPostprocessors]
[vector_pp]
type = ConstantVectorPostprocessor
vector_names = diffusivity_values
value = '1.0 10.0' #we need to set initial values (any values)- these will be over-written
[]
[gradvec]
type = VectorOfPostprocessors
postprocessors = 'grad_bottom grad_top'
execute_on = 'final'
[]
[]

[Executioner]
type = Steady
solve_type = NEWTON
petsc_options_iname = '-pc_type'
petsc_options_value = 'lu'
nl_forced_iters = 1
line_search=none
nl_abs_tol=1e-8
[]

[Outputs]
console = false
file_base = 'adjoint'
[]
Original file line number Diff line number Diff line change
@@ -0,0 +1,139 @@
# This main.i file runs the subapps model.i and grad.i, using an OptimizeFullSolveMultiApp
# The purpose of main.i is to find the two diffusivity_values
# (one in the bottom material of model.i, and one in the top material of model.i)
# such that the misfit between experimental observations (defined in model.i) and MOOSE predictions is minimised.
# The adjoint computed in grad.i is used to compute the gradient for the gradient based LMVM solver in TAO
# PETSc-TAO optimisation is used to perform this inversion
#
[StochasticTools]
[]
[Mesh]
[gmg]
type = GeneratedMeshGenerator
dim = 2
nx = 16
ny = 16
xmin = -4
xmax = 4
ymin = -4
ymax = 4
[]
[bimaterial]
type = SubdomainBoundingBoxGenerator
input = gmg
block_id = 1
bottom_left = '-100 -100 -100'
top_right = '100 0 100'
[]
[name_blocks]
type = RenameBlockGenerator
input = bimaterial
old_block = '0 1'
new_block = 'top bottom'
[]
[]

[OptimizationReporter]
type = ObjectiveGradientMinimize
parameter_names = diffusivity_values
num_values = 2 # diffusivity in the bottom material and in the top material of model.i
initial_condition = '3 4' # the expected result is about '1 10' so this initial condition is not too bad
lower_bounds = '1 1'
upper_bounds = '20 20'
measurement_file = 'synthetic_data.csv'
file_value = 'temperature'
[]

[Executioner]
# type = Optimize
# tao_solver = taoblmvm
# petsc_options_iname = '-tao_fd_gradient -tao_gatol'
# petsc_options_value = ' true 0.001'
type = Optimize
tao_solver = taoblmvm
petsc_options_iname='-tao_gatol'
petsc_options_value='1e-3'
## THESE OPTIONS ARE FOR TESTING THE ADJOINT GRADIENT
# petsc_options_iname='-tao_max_it -tao_fd_test -tao_test_gradient -tao_fd_gradient -tao_fd_delta -tao_gatol'
# petsc_options_value='1 true true false 1e-8 0.1'
# petsc_options = '-tao_test_gradient_view'
# verbose = true
[]
[AuxVariables]
[temperature_forward]
[]
[]
[MultiApps]
[forward]
type = OptimizeFullSolveMultiApp
input_files = model.i
execute_on = "FORWARD"
clone_parent_mesh = true
ignore_solve_not_converge = false
[]
[adjoint]
type = OptimizeFullSolveMultiApp
input_files = grad.i #write this input file to compute the adjoint solution and the gradient
execute_on = "ADJOINT"
clone_parent_mesh = true
ignore_solve_not_converge = false
[]
[]

[Transfers]
[diffusivity_to_forward] #update the model with the new parameters
type = MultiAppReporterTransfer
to_multi_app = forward
from_reporters = 'OptimizationReporter/diffusivity_values'
to_reporters = 'vector_pp/diffusivity_values'
[]
[toForward_measument] #pass the coordinates where we knew the measurements to the forward model to do the extraction of the simulation data at the location of the measurements to compute the misfit
type = MultiAppReporterTransfer
to_multi_app = forward
from_reporters = 'OptimizationReporter/measurement_xcoord OptimizationReporter/measurement_ycoord OptimizationReporter/measurement_zcoord'
to_reporters = 'measure_data/measurement_xcoord measure_data/measurement_ycoord measure_data/measurement_zcoord'
[]
[from_forward] #get the simulation values
type = MultiAppReporterTransfer
from_multi_app = forward
from_reporters = 'data_pt/temperature'
to_reporters = 'OptimizationReporter/simulation_values'
[]
#############
#copy the temprature variable - we will need this for the compuation of the gradient
[fromforwardMesh]
type = MultiAppCopyTransfer
from_multi_app = forward
source_variable = 'temperature'
variable = 'temperature_forward'
[]
[toAdjointMesh]
type = MultiAppCopyTransfer
to_multi_app = adjoint
source_variable = 'temperature_forward'
variable = 'temperature_forward'
[]
#############
[diffusivity_to_adjoint] #update the adjoint with the new parameters
type = MultiAppReporterTransfer
to_multi_app = adjoint
from_reporters = 'OptimizationReporter/diffusivity_values'
to_reporters = 'vector_pp/diffusivity_values'
[]
[toAdjoint]#pas the misfit to the adjoint
type = MultiAppReporterTransfer
to_multi_app = adjoint
from_reporters = 'OptimizationReporter/measurement_xcoord OptimizationReporter/measurement_ycoord OptimizationReporter/measurement_zcoord OptimizationReporter/misfit_values'
to_reporters = 'misfit/measurement_xcoord misfit/measurement_ycoord misfit/measurement_zcoord misfit/misfit_values'
[]
[fromadjoint]
type = MultiAppReporterTransfer
from_multi_app = adjoint
from_reporters = 'gradvec/gradvec'
to_reporters = 'OptimizationReporter/adjoint'
[]
[]

[Outputs]
csv=true
[]
Loading

0 comments on commit 988ae03

Please sign in to comment.