OOFEM  2.4
OOFEM.org - Object Oriented Finite Element Solver
fluidstructureproblem.C
Go to the documentation of this file.
1 /*
2  *
3  * ##### ##### ###### ###### ### ###
4  * ## ## ## ## ## ## ## ### ##
5  * ## ## ## ## #### #### ## # ##
6  * ## ## ## ## ## ## ## ##
7  * ## ## ## ## ## ## ## ##
8  * ##### ##### ## ###### ## ##
9  *
10  *
11  * OOFEM : Object Oriented Finite Element Code
12  *
13  * Copyright (C) 1993 - 2013 Borek Patzak
14  *
15  *
16  *
17  * Czech Technical University, Faculty of Civil Engineering,
18  * Department of Structural Mechanics, 166 29 Prague, Czech Republic
19  *
20  * This library is free software; you can redistribute it and/or
21  * modify it under the terms of the GNU Lesser General Public
22  * License as published by the Free Software Foundation; either
23  * version 2.1 of the License, or (at your option) any later version.
24  *
25  * This program is distributed in the hope that it will be useful,
26  * but WITHOUT ANY WARRANTY; without even the implied warranty of
27  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
28  * Lesser General Public License for more details.
29  *
30  * You should have received a copy of the GNU Lesser General Public
31  * License along with this library; if not, write to the Free Software
32  * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
33  */
34 
35 #include "fluidstructureproblem.h"
36 #include "engngm.h"
37 #include "timestep.h"
38 #include "function.h"
39 #include "metastep.h"
40 #include "exportmodulemanager.h"
41 #include "mathfem.h"
42 #include "oofemtxtdatareader.h"
43 #include "util.h"
44 #include "verbose.h"
45 #include "classfactory.h"
46 
47 #include "pfem.h"
49 #include "../../sm/EngineeringModels/diidynamic.h"
50 #include "../../sm/EngineeringModels/nlineardynamic.h"
51 
52 #include <stdlib.h>
53 
54 #ifdef __OOFEG
55  #include "oofeggraphiccontext.h"
56 #endif
57 
58 namespace oofem {
59 REGISTER_EngngModel(FluidStructureProblem);
60 
62 {
63  ndomains = 1; // domain is needed to store the time step ltf
64 
65  iterationNumber = 0;
66  this->setRenumberFlag();
67 }
68 
70 {}
71 
72 
75 {
76  IRResultType result; // Required by IR_GIVE_FIELD macro
77 
79  if ( result != IRRT_OK ) {
80  return result;
81  }
82 
83  maxiter = 50;
85 
86  rtolv = 1.e-3;
88 
89  rtolp = 1.e-3;
91 
92  return IRRT_OK;
93 }
94 
95 void
97 {
98  for ( auto &emodel: emodelList ) {
99  emodel->initializeYourself(tStep);
100 
101  DIIDynamic *dynamicProblem = dynamic_cast< DIIDynamic * >( emodel.get() );
102  if ( dynamicProblem ) {
104  }
105  NonLinearDynamic *nonlinearDynamicProblem = dynamic_cast< NonLinearDynamic * >( emodel.get() );
106  if ( nonlinearDynamicProblem ) {
107  this->giveCurrentStep()->setTimeDiscretization( nonlinearDynamicProblem->giveInitialTimeDiscretization() );
108  }
109  }
110 }
111 
112 void
114 {
115  PFEM *pfemProblem = dynamic_cast< PFEM * >( this->giveSlaveProblem(1) );
116  Domain *pfemDomain = pfemProblem->giveDomain(1);
117 
118 #ifdef VERBOSE
119  OOFEM_LOG_RELEVANT( "Solving [step number %5d, time %e]\n", stepN->giveNumber(), stepN->giveTargetTime() );
120 #endif
121 
122  FloatArray currentInteractionParticlesVelocities;
123  FloatArray currentInteractionParticlesPressures;
124  FloatArray previousInteractionParticlesVelocities;
125  FloatArray previousInteractionParticlesPressures;
126 
127  if ( stepN->isTheFirstStep() ) {
128  int ndman = pfemDomain->giveNumberOfDofManagers();
129  for ( int i = 1; i <= ndman; i++ ) {
130  if ( dynamic_cast< InteractionPFEMParticle * >( pfemDomain->giveDofManager(i) ) ) {
132  }
133  }
134  }
135 
136  currentInteractionParticlesVelocities.resize( 2 * interactionParticles.giveSize() );
137  previousInteractionParticlesVelocities.resize( 2 * interactionParticles.giveSize() );
138  currentInteractionParticlesPressures.resize( interactionParticles.giveSize() );
139  previousInteractionParticlesPressures.resize( interactionParticles.giveSize() );
140 
141  double velocityDifference = 1.0;
142  double pressureDifference = 1.0;
143 
144  iterationNumber = 0;
145  while ( ( ( velocityDifference > rtolv ) || ( pressureDifference > rtolp ) ) && iterationNumber < maxiter ) {
146  previousInteractionParticlesPressures = currentInteractionParticlesPressures;
147  previousInteractionParticlesVelocities = currentInteractionParticlesVelocities;
148 
149  for ( auto &emodel: emodelList ) {
150  emodel->solveYourselfAt(stepN);
151  }
152 
153  for ( int i = 1; i <= interactionParticles.giveSize(); i++ ) {
154  currentInteractionParticlesPressures.at(i) = pfemProblem->giveUnknownComponent( VM_Total, stepN, pfemDomain, pfemDomain->giveDofManager( interactionParticles.at(i) )->giveDofWithID(P_f) );
155  InteractionPFEMParticle *interactionParticle = dynamic_cast< InteractionPFEMParticle * >( pfemDomain->giveDofManager( interactionParticles.at(i) ) );
156  FloatArray velocities;
157  interactionParticle->giveCoupledVelocities(velocities, stepN);
158  currentInteractionParticlesVelocities.at(2 * ( i - 1 ) + 1) = velocities.at(1);
159  currentInteractionParticlesVelocities.at(2 * ( i - 1 ) + 2) = velocities.at(2);
160  }
161 
162  pressureDifference = 0.0;
163  velocityDifference = 0.0;
164  for ( int i = 1; i <= currentInteractionParticlesPressures.giveSize(); i++ ) {
165  pressureDifference += ( currentInteractionParticlesPressures.at(i) - previousInteractionParticlesPressures.at(i) ) * ( currentInteractionParticlesPressures.at(i) - previousInteractionParticlesPressures.at(i) );
166  velocityDifference += ( currentInteractionParticlesVelocities.at(2 * ( i - 1 ) + 1) - previousInteractionParticlesVelocities.at(2 * ( i - 1 ) + 1) ) * ( currentInteractionParticlesVelocities.at(2 * ( i - 1 ) + 1) - previousInteractionParticlesVelocities.at(2 * ( i - 1 ) + 1) );
167  velocityDifference += ( currentInteractionParticlesVelocities.at(2 * ( i - 1 ) + 2) - previousInteractionParticlesVelocities.at(2 * ( i - 1 ) + 2) ) * ( currentInteractionParticlesVelocities.at(2 * ( i - 1 ) + 2) - previousInteractionParticlesVelocities.at(2 * ( i - 1 ) + 2) );
168  }
169  pressureDifference = sqrt(pressureDifference);
170  velocityDifference = sqrt(velocityDifference);
171 
172  if ( iterationNumber > 0 ) {
173  pressureDifference /= previousInteractionParticlesPressures.computeNorm();
174  velocityDifference /= previousInteractionParticlesVelocities.computeNorm();
175  }
176 
177  OOFEM_LOG_RELEVANT("%3d %le %le\n", iterationNumber++, pressureDifference, velocityDifference);
178  }
179  if ( iterationNumber > maxiter ) {
180  OOFEM_ERROR("Maximal fluid-structure interface iteration count exceded");
181  }
182  stepN->incrementStateCounter();
183 }
184 
185 void
187 {
188  for ( auto &emodel: emodelList ) {
189  emodel->preInitializeNextStep();
190  }
191 }
192 
193 #if 0
194 void
195 FluidStructureProblem :: postInitializeCurrentStep()
196 {
197  for ( int i = 1; i <= nModels; i++ ) {
198  DIIDynamic *dynamicProblem = dynamic_cast< DIIDynamic * >( this->giveSlaveProblem(i) );
199  if ( dynamicProblem ) {
201  }
202  }
203 }
204 #endif
205 } // end namespace oofem
virtual ~FluidStructureProblem()
Destructor.
virtual EngngModel * giveSlaveProblem(int i)
Returns i-th slave problem.
Class and object Domain.
Definition: domain.h:115
int giveNumberOfDofManagers() const
Returns number of dof managers in domain.
Definition: domain.h:432
#define _IFT_FluidStructureProblem_maxiter
double & at(int i)
Coefficient access function.
Definition: floatarray.h:131
#define _IFT_FluidStructureProblem_rtolp
This class represents PFEM method for solving incompressible Navier-Stokes equations.
Definition: pfem.h:124
This class implements nonlinear dynamic engineering problem.
double rtolv
Convergence tolerance.
int iterationNumber
Iteration counter.
double giveTargetTime()
Returns target time.
Definition: timestep.h:146
void incrementStateCounter()
Updates solution state counter.
Definition: timestep.h:190
REGISTER_EngngModel(ProblemSequence)
int & at(int i)
Coefficient access function.
Definition: intarray.h:103
#define OOFEM_LOG_RELEVANT(...)
Definition: logger.h:126
double giveUnknownComponent(ValueModeType mode, TimeStep *tStep, Domain *d, Dof *dof)
Returns requested unknown.
Definition: pfem.C:263
TimeDiscretizationType giveInitialTimeDiscretization()
Definition: diidynamic.h:126
bool isTheFirstStep()
Check if receiver is first step.
Definition: timestep.C:134
void giveCoupledVelocities(FloatArray &answer, TimeStep *stepN)
#define _IFT_FluidStructureProblem_rtolv
int giveNumber()
Returns receiver&#39;s number.
Definition: timestep.h:129
#define OOFEM_ERROR(...)
Definition: error.h:61
int ndomains
Number of receiver domains.
Definition: engngm.h:205
virtual void solveYourselfAt(TimeStep *tStep)
Solves problem for given time step.
Implementation of general sequence (staggered) problem.
virtual IRResultType initializeFrom(InputRecord *ir)
Initializes receiver according to object description in input reader.
FluidStructureProblem(int i, EngngModel *_master=NULL)
Constructor.
Initializes the variable VERBOSE, in order to get a few intermediate messages on screen: beginning an...
Class representing vector of real numbers.
Definition: floatarray.h:82
virtual void preInitializeNextStep()
Does a pre-initialization of the next time step (implement if necessarry)
IRResultType
Type defining the return values of InputRecord reading operations.
Definition: irresulttype.h:47
virtual TimeStep * giveCurrentStep(bool force=false)
Returns current time step.
double computeNorm() const
Computes the norm (or length) of the vector.
Definition: floatarray.C:840
Class representing the general Input Record.
Definition: inputrecord.h:101
This class implements Direct Implicit Integration of Dynamic problem.
Definition: diidynamic.h:73
void followedBy(const IntArray &b, int allocChunk=0)
Appends array b at the end of receiver.
Definition: intarray.C:145
This class represents a fluid particle attached to a node on the structural part of the interface...
TimeDiscretizationType giveInitialTimeDiscretization()
void setTimeDiscretization(TimeDiscretizationType td)
Sets time discretization.
Definition: timestep.h:163
int maxiter
Max number of iterations.
Abstract base class representing the "problem" under consideration.
Definition: engngm.h:181
#define IR_GIVE_OPTIONAL_FIELD(__ir, __value, __id)
Macro facilitating the use of input record reading methods.
Definition: inputrecord.h:78
int giveSize() const
Definition: intarray.h:203
int giveSize() const
Returns the size of receiver.
Definition: floatarray.h:218
the oofem namespace is to define a context or scope in which all oofem names are defined.
virtual void setRenumberFlag()
Sets the renumber flag to true.
Domain * giveDomain(int n)
Service for accessing particular problem domain.
Definition: engngm.C:1720
virtual IRResultType initializeFrom(InputRecord *ir)
Initializes receiver according to object description in input reader.
DofManager * giveDofManager(int n)
Service for accessing particular domain dof manager.
Definition: domain.C:314
std::vector< std::unique_ptr< EngngModel > > emodelList
List of engineering models to solve sequentially.
Class representing solution step.
Definition: timestep.h:80
virtual void initializeYourself(TimeStep *tStep)
Provides the opportunity to initialize state variables stored in element integration points according...
void resize(int s)
Resizes receiver towards requested size.
Definition: floatarray.C:631

This page is part of the OOFEM documentation. Copyright (c) 2011 Borek Patzak
Project e-mail: info@oofem.org
Generated at Tue Jan 2 2018 20:07:28 for OOFEM by doxygen 1.8.11 written by Dimitri van Heesch, © 1997-2011