Skip to content

Commit

Permalink
[SiApp] p-norm damage resp (#12373)
Browse files Browse the repository at this point in the history
* save

* adding p-norm to sensor sum

* test case

* revert files

* apply changes suggested by Suneth

* Update applications/SystemIdentificationApplication/custom_python/add_custom_sensors_to_python.cpp

Co-authored-by: Suneth Warnakulasuriya <[email protected]>

* save

* adding testing

* reverting ref values

---------

Co-authored-by: IAntonau <[email protected]>
Co-authored-by: Suneth Warnakulasuriya <[email protected]>
  • Loading branch information
3 people authored Aug 16, 2024
1 parent 21df993 commit 12643ed
Show file tree
Hide file tree
Showing 12 changed files with 396 additions and 10 deletions.
Original file line number Diff line number Diff line change
Expand Up @@ -99,7 +99,7 @@ void AddCustomSensorsToPython(pybind11::module& m)
;

py::class_<MeasurementResidualResponseFunction, MeasurementResidualResponseFunction::Pointer, AdjointResponseFunction>(sensor_module, "MeasurementResidualResponseFunction")
.def(py::init<>())
.def(py::init<const double>(), py::arg("p_coefficient"))
.def("AddSensor", &MeasurementResidualResponseFunction::AddSensor, py::arg("sensor"))
.def("Clear", &MeasurementResidualResponseFunction::Clear)
.def("GetSensorsList", &MeasurementResidualResponseFunction::GetSensorsList)
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -6,7 +6,8 @@
//
// License: SystemIdentificationApplication/license.txt
//
// Main authors: Suneth Warnakulasuriya
// Main authors: Suneth Warnakulasuriya,
// Ihar Antonau
//

// System includes
Expand Down Expand Up @@ -87,7 +88,8 @@ struct PartialSensitivity

} // namespace MeasurementResidualResponseFunctionUtilities

MeasurementResidualResponseFunction::MeasurementResidualResponseFunction()
MeasurementResidualResponseFunction::MeasurementResidualResponseFunction(const double PCoefficient)
: mPCoefficient(PCoefficient)
{
mResponseGradientList.resize(ParallelUtilities::GetNumThreads());
}
Expand Down Expand Up @@ -147,10 +149,13 @@ double MeasurementResidualResponseFunction::CalculateValue(ModelPart& rModelPart
double value = 0.0;
for (auto& p_sensor : mpSensorsList) {
const double sensor_value = p_sensor->CalculateValue(rModelPart);
value += p_sensor->GetWeight() * std::pow(sensor_value - p_sensor->GetValue(SENSOR_MEASURED_VALUE), 2) * 0.5;
p_sensor->SetSensorValue(sensor_value);
const double current_sensor_error_square = std::pow(sensor_value - p_sensor->GetValue(SENSOR_MEASURED_VALUE), 2) * 0.5;
p_sensor->SetValue(SENSOR_ERROR, current_sensor_error_square);
value += std::pow(p_sensor->GetWeight() * current_sensor_error_square, mPCoefficient);
}
return value;
return std::pow(value, 1 / mPCoefficient);


KRATOS_CATCH("");
}
Expand All @@ -170,11 +175,26 @@ void MeasurementResidualResponseFunction::CalculateDerivative(
rResponseGradient.clear();

auto& local_sensor_response_gradient = mResponseGradientList[OpenMPUtils::ThisThread()];
double temp = 0.0;
for (auto& p_sensor : mpSensorsList) {
temp += ( std::pow( p_sensor->GetValue(SENSOR_ERROR) * 0.5 * p_sensor->GetWeight(), mPCoefficient ) );
}
const double c1 = 1 / mPCoefficient * std::pow( temp, 1/mPCoefficient - 1 );

temp = 0.0;
for (auto& p_sensor : mpSensorsList) {
temp += std::pow( p_sensor->GetWeight() * 0.5 * p_sensor->GetValue(SENSOR_ERROR), mPCoefficient - 1 );
}

const double c2 = mPCoefficient * temp;


for (auto& p_sensor : mpSensorsList) {
TCalculationType::Calculate(*p_sensor, local_sensor_response_gradient, rResidualGradient, rArgs...);
noalias(rResponseGradient) += local_sensor_response_gradient * (p_sensor->GetWeight() * (p_sensor->GetSensorValue() - p_sensor->GetValue(SENSOR_MEASURED_VALUE)));
const double error = std::sqrt( p_sensor->GetValue(SENSOR_ERROR) );
noalias(rResponseGradient) += c1 * c2 * error * local_sensor_response_gradient;
}


KRATOS_CATCH("");
}
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -50,7 +50,7 @@ class KRATOS_API(DIGITAL_TWIN_APPLICATION) MeasurementResidualResponseFunction :
///@{

/// Constructor.
MeasurementResidualResponseFunction();
MeasurementResidualResponseFunction(const double PCoeficient);

/// Destructor.
~MeasurementResidualResponseFunction() override = default;
Expand Down Expand Up @@ -155,6 +155,8 @@ class KRATOS_API(DIGITAL_TWIN_APPLICATION) MeasurementResidualResponseFunction :
///@name Member Variables
///@{

double mPCoefficient;

std::vector<Sensor::Pointer> mpSensorsList;

std::vector<Vector> mResponseGradientList;
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -199,4 +199,4 @@ def _UpdateAndOutputFields(self, update: ContainerExpressionTypes) -> None:
un_buffered_data.SetValue(f"{self.controlled_physical_variable.Name()}_physical_phi_derivative", self.physical_phi_derivative_field.Clone(), overwrite=True)

def __str__(self) -> str:
return f"Control [type = {self.__class__.__name__}, name = {self.GetName()}, model part name = {self.adjoint_model_part_operation.GetModelPartFullName()}, control variable = {self.controlled_physical_variable.Name()}"
return f"Control [type = {self.__class__.__name__}, name = {self.GetName()}, model part name = {self.adjoint_model_part_operation.GetModelPartFullName()}, control variable = {self.controlled_physical_variable.Name()}"
Original file line number Diff line number Diff line change
Expand Up @@ -92,7 +92,7 @@ def Initialize(self) -> None:

def Check(self) -> None:
pass

def Finalize(self) -> None:
self.adjoint_analysis.Finalize()

Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -22,6 +22,7 @@ def Initialize(self):
"perturbation_size" : 1e-8,
"adapt_perturbation_size" : true,
"list_of_sensors" : [],
"p_coefficient" : 1,
"output_settings" : {
"output_sensor_sensitivity_fields": false,
"output_folder" : "Optimization_Results/sensor_sensitivity_fields"
Expand All @@ -36,7 +37,9 @@ def Initialize(self):
model_part.ProcessInfo[KratosSI.ADAPT_PERTURBATION_SIZE] = sensor_settings["adapt_perturbation_size"].GetBool()
self.listof_sensors = GetSensors(model_part, sensor_settings["list_of_sensors"].values())

self.measurement_residual_response_function = KratosSI.Sensors.MeasurementResidualResponseFunction()
p_coefficient = sensor_settings["p_coefficient"].GetDouble()
self.measurement_residual_response_function = KratosSI.Sensors.MeasurementResidualResponseFunction(p_coefficient)

for sensor in self.listof_sensors:
sensor.SetValue(KratosSI.SENSOR_MEASURED_VALUE, 0.0)
self.measurement_residual_response_function.AddSensor(sensor)
Expand Down
Original file line number Diff line number Diff line change
@@ -0,0 +1,3 @@
{
"properties": []
}
Original file line number Diff line number Diff line change
@@ -0,0 +1,77 @@
{
"problem_data": {
"problem_name": "benchmark_2",
"parallel_type": "OpenMP",
"start_time": 0,
"end_time": 1,
"echo_level": 0
},
"sensor_settings": {
"perturbation_size": 1e-8,
"adapt_perturbation_size": true,
"p_coefficient": 4,
"@include_json": "auxiliary_files/system_identification/sensor_data.json"
},
"solver_settings": {
"solver_type": "adjoint_static",
"analysis_type": "linear",
"model_part_name": "AdjointStructure",
"domain_size": 3,
"time_stepping": {
"time_step": 1.0
},
"compute_reactions": false,
"move_mesh_flag": false,
"sensitivity_settings": {
"sensitivity_model_part_name": "all_nodes_elements_model_part",
"element_data_value_sensitivity_variables": [
"YOUNG_MODULUS"
],
"build_mode": "static"
},
"echo_level": 0,
"rotation_dofs": true,
"model_import_settings": {
"input_type": "use_input_model_part"
},
"material_import_settings": {
"materials_filename": "auxiliary_files/system_identification/adjoint_material_parameters.json"
},
"linear_solver_settings": {
"solver_type": "amgcl"
}
},
"processes": {
"constraints_process_list": [
{
"python_module": "assign_vector_variable_process",
"kratos_module": "KratosMultiphysics",
"help": "This process fixes the selected components of a given vector variable",
"process_name": "AssignVectorVariableProcess",
"Parameters": {
"mesh_id": 0,
"model_part_name": "AdjointStructure.BoundPtsType7",
"variable_name": "ADJOINT_DISPLACEMENT",
"value": [
0.0,
0.0,
0.0
],
"constrained": [
true,
true,
true
],
"interval": [
0.0,
"End"
]
}
}
],
"loads_process_list": [],
"list_other_processes": []
},
"output_processes": {},
"analysis_stage": "KratosMultiphysics.SystemIdentificationApplication.sensor_sensitivity_solvers.sensor_sensitivity_static_analysis"
}
Original file line number Diff line number Diff line change
@@ -0,0 +1,133 @@
{
"problem_data": {
"parallel_type": "OpenMP",
"echo_level": 0
},
"model_parts": [
{
"type": "mdpa_model_part_controller",
"settings": {
"model_part_name": "AdjointStructure",
"domain_size": 3,
"input_filename": "auxiliary_files/structure"
}
},
{
"type": "connectivity_preserving_model_part_controller",
"settings": {
"transformation_settings": [
{
"source_model_part_name": "AdjointStructure",
"destination_model_part_name": "Structure",
"destination_element_name": "ShellThinElement3D3N",
"destination_condition_name": "LineLoadCondition3D2N"
}
]
}
}
],
"analyses": [
{
"name": "Structure_static",
"type": "kratos_analysis_execution_policy",
"settings": {
"model_part_names": [
"Structure"
],
"analysis_module": "KratosMultiphysics.StructuralMechanicsApplication",
"analysis_type": "StructuralMechanicsAnalysis",
"analysis_settings": {
"@include_json": "auxiliary_files/system_identification/primal_project_parameters.json"
}
}
}
],
"responses": [
{
"name": "damage_response",
"type": "damage_detection_response",
"module": "KratosMultiphysics.SystemIdentificationApplication.responses",
"settings": {
"evaluated_model_part_names": [
"AdjointStructure"
],
"adjoint_parameters": {
"@include_json": "auxiliary_files/system_identification_p_norm/adjoint_project_parameters.json"
},
"test_analysis_list": [
{
"primal_analysis_name": "Structure_static",
"sensor_measurement_csv_file": "auxiliary_files/damaged_problem/measured_data_ref.csv",
"weight": 1.0
}
]
}
}
],
"controls": [
{
"name": "material_control",
"type": "material_properties_control",
"module": "KratosMultiphysics.SystemIdentificationApplication.controls",
"settings": {
"model_part_names": [
{
"primal_model_part_name": "Structure",
"adjoint_model_part_name": "AdjointStructure"
}
],
"control_variable_name": "YOUNG_MODULUS",
"control_variable_bounds": [
0.0,
30000000000.0
],
"filter_settings": {
"filter_type": "explicit_filter",
"filter_radius_settings": {
"filter_radius": 5.0,
"filter_radius_type": "constant"
}
}
}
}
],
"algorithm_settings": {
"type": "algorithm_steepest_descent",
"settings": {
"echo_level": 0,
"line_search": {
"type": "const_step",
"init_step": 0.01,
"gradient_scaling": "inf_norm"
},
"conv_settings": {
"type": "max_iter",
"max_iter": 5
}
},
"controls": [
"material_control"
],
"objective": {
"response_name": "damage_response",
"type": "minimization",
"scaling": 1.0
}
},
"processes": {
"kratos_processes": {},
"optimization_data_processes": {
"output_processes": [
{
"type": "optimization_problem_ascii_output_process",
"module": "KratosMultiphysics.OptimizationApplication.processes",
"settings": {
"output_file_name": "auxiliary_files/summary_p_norm.csv",
"write_kratos_version": false,
"write_time_stamp": false
}
}
]
}
}
}
Loading

0 comments on commit 12643ed

Please sign in to comment.