Kratos: [Core][Discussion] Adaptive strategy

Created on 5 Apr 2019  路  24Comments  路  Source: KratosMultiphysics/Kratos

i think it would be useful to implement an adaptive time stepping in the analysis stages of the structure (although it may be useful also for the fluid)

a tentative implementation could be something on the line of the following code (the implementation was done in a rush and is not good quality...but the idea is there):

~~~py
from __future__ import print_function, absolute_import, division #makes KratosMultiphysics backward compatible with python 2.6 and 2.7

import KratosMultiphysics
from KratosMultiphysics.StructuralMechanicsApplication.structural_mechanics_analysis import StructuralMechanicsAnalysis

"""
For user-scripting it is intended that a new class is derived
from StructuralMechanicsAnalysis to do modifications
"""

class ModifiedStructuralMechanicsAnalysis(StructuralMechanicsAnalysis):

def __init__(self,model,project_parameters):
    super(StructuralMechanicsAnalysis,self).__init__(model,project_parameters)

    self.reduction_factor = 0.5
    self.incrementation_factor = 1.1
    self.max_delta_time = project_parameters["solver_settings"]["time_stepping"]["time_step"].GetDouble()
    self.delta_time = self.max_delta_time
    self.max_allowed_iterations = 10

def FinalizeSolutionStep(self):
    super(StructuralMechanicsAnalysis,self).FinalizeSolutionStep()

    max_allowed_iterations = 10
    if(self._GetSolver().GetComputingModelPart().ProcessInfo[KratosMultiphysics.NL_ITERATION_NUMBER] >= max_allowed_iterations):
        raise Exception("max_number_of_iterations_exceeded")

def RunSolutionLoop(self):
    """This function executes the solution loop of the AnalysisStage
    It can be overridden by derived classes
    """
    while self.KeepAdvancingSolutionLoop():
        self.delta_time *= self.incrementation_factor
        if(self.delta_time > self.max_delta_time):
            self.delta_time = self.max_delta_time
        t = self._GetSolver().GetComputingModelPart().ProcessInfo[KratosMultiphysics.TIME]
        new_time = t + self.delta_time
        self._GetSolver().GetComputingModelPart().ProcessInfo[KratosMultiphysics.STEP] += 1
        self._GetSolver().main_model_part.CloneTimeStep(new_time)

        converged = False
        while(not converged):

            t = self._GetSolver().GetComputingModelPart().ProcessInfo[KratosMultiphysics.TIME]
            corrected_time = t - self._GetSolver().GetComputingModelPart().ProcessInfo[KratosMultiphysics.DELTA_TIME] + self.delta_time
            self._GetSolver().GetComputingModelPart().ProcessInfo[KratosMultiphysics.TIME] = corrected_time
            self._GetSolver().GetComputingModelPart().ProcessInfo[KratosMultiphysics.DELTA_TIME] = self.delta_time

            print("----------------------- t = -------------------------",self._GetSolver().GetComputingModelPart().ProcessInfo[KratosMultiphysics.TIME])

            self.InitializeSolutionStep()
            self._GetSolver().Predict()
            self._GetSolver().SolveSolutionStep()

            if(self._GetSolver().GetComputingModelPart().ProcessInfo[KratosMultiphysics.NL_ITERATION_NUMBER] >= self.max_allowed_iterations):
                self.delta_time *= self.reduction_factor
                converged = False
                for node in self._GetSolver().GetComputingModelPart().Nodes:
                    dold = node.GetSolutionStepValue(KratosMultiphysics.DISPLACEMENT,1)
                    node.SetSolutionStepValue(KratosMultiphysics.DISPLACEMENT,0,dold)
                    node.X = node.X0 + dold[0]
                    node.Y = node.Y0 + dold[1]
                    node.Z = node.Z0 + dold[2]
                print("***************************************************+")
                print("***************************************************+")
                print("***************************************************+")
                print("***************************************************+")
            else:
                converged = True
                self.delta_time *= self.incrementation_factor

        self.FinalizeSolutionStep()
        self.OutputSolutionStep()

if __name__ == "__main__":

with open("ProjectParameters.json",'r') as parameter_file:
    parameters = KratosMultiphysics.Parameters(parameter_file.read())

model = KratosMultiphysics.Model()
simulation = ModifiedStructuralMechanicsAnalysis(model,parameters)
simulation.Run()

~~~

Discussion

All 24 comments

Nice :+1:
This will be very handy

This can be simplified very much after #4535

BTW we should be really careful with what to integrate in the AnalysisStages of the apps
E.g. this is highly conflictive with #4558

Also note that we should avoid using the computingmodelpart as much as possible, since

  1. It will most likely eventually disappear
  2. It does not make sense for coupled solvers

I don't like the prints...

Apart from that I agrees with @philbucher

I don't like the prints...

I think this is just a dummy implementation as a proof-of-concept, that won't go into master :)

I would suggest to use the ModelPart::ReduceTimeStep instead of manually reset the displacement. Meanwhile we should have a process which reverts the nodal coordinates.

The problem with ReduceTimeStep (which i agree using) is that it does not reset the X position

so one needs to manually do X = X0 + DISPLACEMENT_X

anyway, anyone willing to take the lead on this? i really cannot right now

The problem with ReduceTimeStep (which i agree using) is that it does not reset the X position

so one needs to manually do X = X0 + DISPLACEMENT_X

@rubenzorrilla recently created an utility on VariableUtils for this

The problem with ReduceTimeStep (which i agree using) is that it does not reset the X position
so one needs to manually do X = X0 + DISPLACEMENT_X

@rubenzorrilla recently created an utility on VariableUtils for this

Yes, and can think on enhance it (or do something similar) to do the operation @RiccardoRossi is pointing out.

What is the status of this? Is it already implemented in Kratos?

I am waiting for the strategies factory PR to be merged

Can we move the part, which resets the displacements to C++? I guess it would be time-consuming for large models to be done in python.

Can we move the part, which resets the displacements to C++? I guess it would be time-consuming for large models to be done in python.

It is already in the VariableUtils

Can we move the part, which resets the displacements to C++? I guess it would be time-consuming for large models to be done in python.

It is already in the VariableUtils

If I'm not wrong, what you ask for is the UpdateCurrentPosition function of the VariableUtils.

I invoke @KratosMultiphysics/implementation-committee to discuss this further

Thanks!

I checked all the functions in VariableUtils but I think none of them is doing the job. Maybe I am missing something but what is happining in UpdateCurrentPosition is X = X0 + DISPLACEMENT_X_current while here we want to have X = X0 + DISPLACEMENT_X_old

One sec

@Vahid-Galavi I updated the function in #5432

The @KratosMultiphysics/implementation-committee prompts @RiccardoRossi to have a look at #4558, which is, as mentioned by @philbucher in conflict with this proposal.

ping @RiccardoRossi

ping @RiccardoRossi

And remember to review first #4558

On behalf of @KratosMultiphysics/implementation-committee

Considering the usefulness of this, we think that an abstract solution could be reached and implemented in the base analysis stage. The new time step can be obtained by a function like AdjustTimeStep or AdaptTimeStep and the resetting of dofs can be done via Scheme in C++ Eg: through a flag or a setting. One can implement it in an Application as a prototype and then its easier everyone to review it.

On behalf of @KratosMultiphysics/implementation-committee

Considering the usefulness of this, we think that an abstract solution could be reached and implemented in the base analysis stage. The new time step can be obtained by a function like AdjustTimeStep or AdaptTimeStep and the resetting of dofs can be done via Scheme in C++ Eg: through a flag or a setting. One can implement it in an Application as a prototype and then its easier everyone to review it.

Remember to review first #4558, which conflicts with this

Was this page helpful?
0 / 5 - 0 ratings

Related issues

armingeiser picture armingeiser  路  6Comments

riccardotosi picture riccardotosi  路  5Comments

Vahid-Galavi picture Vahid-Galavi  路  4Comments

josep-m-carbonell picture josep-m-carbonell  路  4Comments

roigcarlo picture roigcarlo  路  6Comments