LifeV
OseenSolverShapeDerivative.hpp
Go to the documentation of this file.
1 //@HEADER
2 /*
3 *******************************************************************************
4 
5  Copyright (C) 2004, 2005, 2007 EPFL, Politecnico di Milano, INRIA
6  Copyright (C) 2010 EPFL, Politecnico di Milano, Emory University
7 
8  This file is part of LifeV.
9 
10  LifeV is free software; you can redistribute it and/or modify
11  it under the terms of the GNU Lesser General Public License as published by
12  the Free Software Foundation, either version 3 of the License, or
13  (at your option) any later version.
14 
15  LifeV is distributed in the hope that it will be useful,
16  but WITHOUT ANY WARRANTY; without even the implied warranty of
17  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
18  Lesser General Public License for more details.
19 
20  You should have received a copy of the GNU Lesser General Public License
21  along with LifeV. If not, see <http://www.gnu.org/licenses/>.
22 
23 *******************************************************************************
24 */
25 //@HEADER
26 
27 /*!
28  @file
29  @brief This file contains an Oseen equation solver class with shape derivative
30  for fluid structure interaction problem
31 
32  @author MiGlobaluel A. Fernandez <miGlobaluel.fernandez@inria.fr>
33  Christoph Winkelmann <christoph.winkelmann@epfl.ch>
34  @date 09-06-2003
35 
36  @author G. Fourestey
37  @date 00-02-2007
38 
39  @contributor Zhen Wang <zhen.wang@emory.edu>
40 
41  */
42 
43 
44 #ifndef OSEENSOLVERSHAPEDERIVATIVE_H
45 #define OSEENSOLVERSHAPEDERIVATIVE_H 1
46 
47 #include <lifev/navier_stokes/solver/OseenSolver.hpp>
48 
49 namespace LifeV
50 {
51 
52 //! @class OseenSolverShapeDerivative
53 /*!
54 
55  @brief This class contains an Oseen equation solver class with shape derivative
56  for fluid structure interaction problem
57 
58  @author MiGlobaluel A. Fernandez <miGlobaluel.fernandez@inria.fr>
59  Christoph Winkelmann <christoph.winkelmann@epfl.ch>
60  @date 09-06-2003
61  @author G. Fourestey
62  @date 02-2007
63 
64  @contributor Zhen Wang <zhen.wang@emory.edu>
65 
66  */
67 
68 template< typename MeshType, typename SolverType = LifeV::SolverAztecOO >
69 class OseenSolverShapeDerivative:
70  public OseenSolver< MeshType, SolverType >
71 {
72 
73 public:
74 
75  //! @name Public Types
76  //@{
77 
78  typedef MeshType mesh_Type;
79  typedef SolverType linearSolver_Type;
80  typedef OseenSolver< mesh_Type, linearSolver_Type > oseenSolver_Type;
81  typedef typename oseenSolver_Type::vector_Type vector_Type;
82  typedef typename oseenSolver_Type::solution_Type solution_Type;
83  typedef typename oseenSolver_Type::solutionPtr_Type solutionPtr_Type;
84  typedef typename oseenSolver_Type::matrix_Type matrix_Type;
85  typedef typename oseenSolver_Type::matrixPtr_Type matrixPtr_Type;
86  typedef typename oseenSolver_Type::data_Type data_Type;
87  typedef typename oseenSolver_Type::preconditioner_Type preconditioner_Type;
88  typedef typename oseenSolver_Type::preconditionerPtr_Type preconditionerPtr_type;
89  typedef typename oseenSolver_Type::bcHandler_Type bcHandler_Type;
90 
91  //@}
92 
93  //! @name Constructors & Destructor
94  //@{
95 
96  //! Empty constructor
97  OseenSolverShapeDerivative();
98 
99  //! Constructor
100  /*!
101  @param dataType OseenData class
102  @param velocityFESpace Velocity FE space
103  @param pressureFESpace Pressure FE space
104  @param communicator MPI communicator
105  @param lagrangeMultiplier Lagrange multiplier
106  */
107  OseenSolverShapeDerivative ( std::shared_ptr<data_Type> dataType,
108  FESpace<mesh_Type, MapEpetra>& velocityFESpace,
109  FESpace<mesh_Type, MapEpetra>& pressureFESpace,
110  std::shared_ptr<Epetra_Comm>& communicator,
111  const Int lagrangeMultiplier = 0);
112 
113  //! Constructor
114  /*!
115  @param dataType OseenData class
116  @param velocityFESpace Velocity FE space
117  @param pressureFESpace Pressure FE space
118  @param communicator MPI communicator
119  @param monolithicMap MapEpetra class
120  @param offset
121  */
122  OseenSolverShapeDerivative ( std::shared_ptr<data_Type> dataType,
123  FESpace<mesh_Type, MapEpetra>& velocityFESpace,
124  FESpace<mesh_Type, MapEpetra>& pressureFESpace,
125  std::shared_ptr<Epetra_Comm>& communicator,
126  const MapEpetra monolithicMap,
127  const UInt offset = 0 );
128 
129  //! Constructor
130  /*!
131  @param dataType OseenData class
132  @param velocityFESpace Velocity FE space
133  @param pressureFESpace Pressure FE space
134  @param mmFESpace FE space
135  @param communicator MPI communicator
136  @param monolithicMap MapEpetra class
137  @param offset
138  */
139  OseenSolverShapeDerivative ( std::shared_ptr<data_Type> dataType,
140  FESpace<mesh_Type, MapEpetra>& velocityFESpace,
141  FESpace<mesh_Type, MapEpetra>& pressureFESpace,
142  FESpace<mesh_Type, MapEpetra>& mmFESpace,
143  std::shared_ptr<Epetra_Comm>& communicator,
144  const MapEpetra monolithicMap,
145  const UInt offset = 0 );
146 
147  //! Virtual destructor
148  virtual ~OseenSolverShapeDerivative();
149 
150  //@}
151 
152 
153  //! @name Methods
154  //@{
155 
156  //! Set up data from GetPot
157  /*!
158  @param dataFile GetPot object
159  */
160  void setUp ( const GetPot& dataFile );
161 
162  //!
163  /*!
164  @param bcHandler BC handler
165  */
166  void solveLinearSystem ( bcHandler_Type& bcHandler );
167 
168  //! Update linear system.
169  /*!
170  @param matrixNoBC Fluid matrix withoud BC
171  @param alpha alpha
172  @param un Beta
173  @param uk Fluid solution
174  @param disp mesh_Type deltaX
175  @param w mesh_Type Velocity
176  @param dw mesh_Type deltaVelocity
177  @param sourceVector RHS (usually 0 )
178  */
179  void updateLinearSystem ( const matrix_Type& matrixNoBC,
180  Real& alpha,
181  const vector_Type& un,
182  const vector_Type& uk,
183  const vector_Type& disp,
184  const vector_Type& w,
185  const vector_Type& dw,
186  const vector_Type& sourceVector);
187 
188  //! Update shape derivatives.
189  /*!
190  @param matrixNoBC Fluid matrix withoud BC
191  @param alpha alpha
192  @param un Beta
193  @param uk Fluid solution
194  @param disp mesh_Type deltaX
195  @param w mesh_Type Velocity
196  @param offset
197  @param dFESpace
198  @param wImplicit
199  @param convectiveTermDerivative
200  */
201  void updateShapeDerivatives ( matrix_Type& matrixNoBC,
202  Real& alpha,
203  const vector_Type& un,
204  const vector_Type& uk,
205  //const vector_Type& disp,
206  const vector_Type& w,
207  UInt offset,
208  FESpace<mesh_Type, MapEpetra>& dFESpace,
209  bool wImplicit = true,
210  bool convectiveTermDerivative = false);
211 
212  //@}
213 
214 
215  //! @name Set Methods
216  //@{
217 
218  //! Set
219  /*!
220  @param rightHandSide
221  */
222  void updateLinearRightHandSideNoBC ( const vector_Type& rightHandSide)
223  {
224  M_linearRightHandSideNoBC = rightHandSide;
225  }
226 
227  //@}
228 
229  //! @name Get Methods
230  //@{
231 
232  //! Return
233  /*!
234  @return M_linearRightHandSideNoBC
235  */
236  vector_Type& linearRightHandSideNoBC()
237  {
238  return M_linearRightHandSideNoBC;
239  }
240 
241  const vector_Type& linearRightHandSideNoBC() const
242  {
243  return M_linearRightHandSideNoBC;
244  }
245 
246  //! Get the solution of the Shape Derivative problem.
247  /*!
248  @return vector containing the solution of the Shape Derivative problem.
249  */
250  const vector_Type& linearSolution() const
251  {
252  return M_linearSolution;
253  }
254 
255  //! Return
256  /*!
257  @return stabilization
258  */
259  bool stabilization()
260  {
261  return this->M_stabilization;
262  }
263 
264  const bool& stabilization() const
265  {
266  return this->M_stabilization;
267  }
268 
269  //! Compute the derivative of the flow rate on a boundary face
270  /*!
271  * @param flag boundary flag
272  * @return derivative of flow rate
273  */
274  Real linearFlux ( const markerID_Type& flag );
275  LIFEV_DEPRECATED ( Real getLinearFlux ( const markerID_Type& flag ) );
276 
277  //! Compute the derivative of the pressure on a boundary face
278  /*!
279  * @param flag boundary flag
280  * @return derivative of pressure
281  */
282  Real linearPressure ( const markerID_Type& flag );
283  LIFEV_DEPRECATED ( Real getLinearPressure ( const markerID_Type& flag ) );
284 
285  //! Compute the derivative of a Lagrange multiplier (which correspond to the value of the derivative of the mean normal stress on a boundary face)
286  /*!
287  * @param flag flag of the boundary face associated with the flux and the Lagrange multiplier we want.
288  * @param BC BChandler containing the boundary conditions of the problem.
289  * @return derivative of the Lagrange multiplier
290  */
291  Real linearLagrangeMultiplier ( const markerID_Type& flag, bcHandler_Type& bcHandler );
292  LIFEV_DEPRECATED ( Real getLinearLagrangeMultiplier ( const markerID_Type& flag, bcHandler_Type& bcHandler ) );
293 
294  //! Compute the derivative of the kinetic normal stress (i.e., the normal stress due to the kinetic energy) on a boundary face with given flag
295  /*!
296  * @see \cite BlancoMalossi2012 \cite Malossi-Thesis
297  *
298  * @param flag boundary flag
299  * @return derivative of the kinetic normal stress
300  */
301  Real linearKineticNormalStress ( const markerID_Type& flag );
302 
303  //! Compute the derivative of the kinetic normal stress (i.e., the normal stress due to the kinetic energy) on a boundary face with a given flag and a given solution
304  /*!
305  * @see \cite BlancoMalossi2012 \cite Malossi-Thesis
306  *
307  * @param flag boundary flag
308  * @param solution problem solution
309  * @param linearSolution linear problem solution
310  * @return derivative of the kinetic normal stress
311  */
312  Real linearKineticNormalStress ( const markerID_Type& flag, const vector_Type& solution, const vector_Type& linearSolution );
313 
314  //! Compute the derivative of the mean normal stress on a boundary face with a given flag
315  /*!
316  * @see \cite BlancoMalossi2012 \cite Malossi-Thesis
317  *
318  * @param flag flag of the boundary face associated with the flux and the Lagrange multiplier we want.
319  * @param BC BChandler containing the boundary conditions of the problem.
320  * @return derivative of the mean normal stress
321  */
322  Real linearMeanNormalStress ( const markerID_Type& flag, bcHandler_Type& bcHandler );
323 
324  //! Compute the derivative of the mean normal stress on a boundary face with a given flag
325  /*!
326  * @see \cite BlancoMalossi2012 \cite Malossi-Thesis
327  *
328  * @param flag flag of the boundary face associated with the flux and the Lagrange multiplier we want.
329  * @param BC BChandler containing the boundary conditions of the problem.
330  * @param linearSolution linear problem solution
331  * @return derivative of the mean normal stress
332  */
333  Real linearMeanNormalStress ( const markerID_Type& flag, bcHandler_Type& bcHandler, const vector_Type& linearSolution );
334 
335  //! Compute the derivative of the mean total normal stress on a boundary face with a given flag
336  /*!
337  * @see \cite BlancoMalossi2012 \cite Malossi-Thesis
338  *
339  * @param flag flag of the boundary face associated with the flux and the Lagrange multiplier we want.
340  * @param BC BChandler containing the boundary conditions of the problem.
341  * @return derivative of the mean total normal stress
342  */
343  Real linearMeanTotalNormalStress ( const markerID_Type& flag, bcHandler_Type& bcHandler );
344 
345  //! Compute the derivative of the mean total normal stress on a boundary face with a given flag
346  /*!
347  * @see \cite BlancoMalossi2012 \cite Malossi-Thesis
348  *
349  * @param flag flag of the boundary face associated with the flux and the Lagrange multiplier we want.
350  * @param BC BChandler containing the boundary conditions of the problem.
351  * @param solution problem solution
352  * @param linearSolution linear problem solution
353  * @return derivative of the mean total normal stress
354  */
355  Real linearMeanTotalNormalStress ( const markerID_Type& flag, bcHandler_Type& bcHandler, const vector_Type& solution, const vector_Type& linearSolution );
356 
357  //@}
358 
359 
360 private:
361 
362  //@{
363 
364  //! Empty copy constructor
365  OseenSolverShapeDerivative ( const OseenSolverShapeDerivative& oseenShapeDerivative );
366 
367  //@}
368 
369  vector_Type M_linearRightHandSideNoBC;
370  vector_Type M_linearRightHandSideFull;
371 
372  vector_Type M_linearSolution;
373 
374  linearSolver_Type M_linearLinSolver;
375  preconditionerPtr_type M_linearPreconditioner;
376 
377 
378  VectorElemental M_elementVectorVelocity; // Elementary right hand side for the linearized velocity
379  VectorElemental M_elementVectorPressure; // Elementary right hand side for the linearized pressure
380  // std::shared_ptr<MatrixElemental> M_elementMatrixVelocity;
381  // std::shared_ptr<MatrixElemental> M_elementMatrixConvective;
382  // std::shared_ptr<MatrixElemental> M_elementMatrixPressure; // Elementary displacement for right hand side
383  VectorElemental M_elementMeshVelocity; // Elementary mesh velocity
384  VectorElemental M_elementVelocity; // Elementary velocity
385  VectorElemental M_elementPressure; // Elementary pressure
386  VectorElemental M_elementConvectionVelocity; // Elementary convection velocity
387  VectorElemental M_elementDisplacement; // Elementary displacement for right hand side
388  VectorElemental M_elementVelocityRightHandSide; // Elementary mesh velocity for right hand side
389  VectorElemental M_u_loc;
390  bool M_reuseLinearPreconditioner;
391  FESpace<mesh_Type, MapEpetra>* M_mmFESpace;
392 };
393 
394 
395 
396 // ===================================================
397 // Constructors & Destructor
398 // ===================================================
399 
400 template<typename MeshType, typename SolverType>
401 OseenSolverShapeDerivative<MeshType, SolverType>::
402 OseenSolverShapeDerivative ( std::shared_ptr<data_Type> dataType,
403  FESpace<mesh_Type, MapEpetra>& velocityFESpace,
404  FESpace<mesh_Type, MapEpetra>& pressureFESpace,
405  std::shared_ptr<Epetra_Comm>& communicator,
406  const Int lagrangeMultiplier) :
407  oseenSolver_Type ( dataType,
408  velocityFESpace,
409  pressureFESpace,
410  communicator,
411  lagrangeMultiplier ),
412  M_linearRightHandSideNoBC ( this->getMap() ),
413  M_linearRightHandSideFull ( this->getMap() ),
414  M_linearSolution ( this->getMap() ),
415  M_linearLinSolver ( communicator ),
416  M_linearPreconditioner ( ),
417  M_elementVectorVelocity ( this->M_velocityFESpace.fe().nbFEDof(), nDimensions ),
418  M_elementVectorPressure ( this->M_pressureFESpace.fe().nbFEDof(), 1 ),
419  // M_elementVectorPressure ( this->M_pressureFESpace.fe().nbFEDof(), nDimensions ),
420  M_elementMeshVelocity ( this->M_velocityFESpace.fe().nbFEDof(), nDimensions ),
421  M_elementVelocity ( this->M_velocityFESpace.fe().nbFEDof(), nDimensions ),
422  M_elementPressure ( this->M_pressureFESpace.fe().nbFEDof(), 1 ),
423  M_elementConvectionVelocity ( this->M_velocityFESpace.fe().nbFEDof(), nDimensions ),
424  M_elementDisplacement ( this->M_velocityFESpace.fe().nbFEDof(), nDimensions ),
425  M_elementVelocityRightHandSide ( this->M_velocityFESpace.fe().nbFEDof(), nDimensions ),
426  M_u_loc ( this->M_velocityFESpace.fe().nbFEDof(), nDimensions ),
427  M_reuseLinearPreconditioner ( true )
428 {
429 
430 }
431 
432 template<typename MeshType, typename SolverType>
433 OseenSolverShapeDerivative<MeshType, SolverType>::
434 OseenSolverShapeDerivative ( std::shared_ptr<data_Type> dataType,
435  FESpace<mesh_Type, MapEpetra>& velocityFESpace,
436  FESpace<mesh_Type, MapEpetra>& pressureFESpace,
437  std::shared_ptr<Epetra_Comm>& communicator,
438  const MapEpetra monolithicMap,
439  const UInt offset ) :
440  oseenSolver_Type ( dataType,
441  velocityFESpace,
442  pressureFESpace,
443  communicator,
444  monolithicMap,
445  offset),
446  M_linearRightHandSideNoBC ( this->getMap() ),
447  M_linearRightHandSideFull ( this->getMap() ),
448  M_linearSolution ( this->getMap() ),
449  M_linearLinSolver ( communicator ),
450  M_linearPreconditioner ( ),
451  M_elementVectorVelocity ( this->M_velocityFESpace.fe().nbFEDof(), nDimensions ),
452  M_elementVectorPressure ( this->M_pressureFESpace.fe().nbFEDof(), 1 ),
453  // M_elementVectorPressure ( this->M_pressureFESpace.fe().nbFEDof(), nDimensions ),
454  M_elementMeshVelocity ( this->M_velocityFESpace.fe().nbFEDof(), nDimensions ),
455  M_elementVelocity ( this->M_velocityFESpace.fe().nbFEDof(), nDimensions ),
456  M_elementPressure ( this->M_pressureFESpace.fe().nbFEDof(), 1 ),
457  M_elementConvectionVelocity ( this->M_velocityFESpace.fe().nbFEDof(), nDimensions ),
458  M_elementDisplacement ( this->M_velocityFESpace.fe().nbFEDof(), nDimensions ),
459  M_elementVelocityRightHandSide ( this->M_velocityFESpace.fe().nbFEDof(), nDimensions ),
460  M_u_loc ( this->M_velocityFESpace.fe().nbFEDof(), nDimensions ),
461  M_reuseLinearPreconditioner ( true )
462 {
463 
464 }
465 
466 template<typename MeshType, typename SolverType>
467 OseenSolverShapeDerivative<MeshType, SolverType>::
468 OseenSolverShapeDerivative ( std::shared_ptr<data_Type> dataType,
469  FESpace<mesh_Type, MapEpetra>& velocityFESpace,
470  FESpace<mesh_Type, MapEpetra>& pressureFESpace,
471  FESpace<mesh_Type, MapEpetra>& mmFESpace,
472  std::shared_ptr<Epetra_Comm>& communicator,
473  const MapEpetra monolithicMap,
474  const UInt offset) :
475  oseenSolver_Type (dataType,
476  velocityFESpace,
477  pressureFESpace,
478  communicator,
479  monolithicMap,
480  offset),
481  M_linearRightHandSideNoBC ( this->getMap() ),
482  M_linearRightHandSideFull ( this->getMap() ),
483  M_linearSolution ( this->getMap() ),
484  M_linearLinSolver( ),
485  M_linearPreconditioner ( ),
486  M_elementVectorVelocity ( this->M_velocityFESpace.fe().nbFEDof(), nDimensions ),
487  M_elementVectorPressure ( this->M_pressureFESpace.fe().nbFEDof(), 1 ),
488  // M_elementVectorPressure ( this->M_pressureFESpace.fe().nbFEDof(), nDimensions ),
489  M_elementMeshVelocity ( this->M_velocityFESpace.fe().nbFEDof(), nDimensions ),
490  M_elementVelocity ( this->M_velocityFESpace.fe().nbFEDof(), nDimensions ),
491  M_elementPressure ( this->M_pressureFESpace.fe().nbFEDof(), 1 ),
492  M_elementConvectionVelocity ( this->M_velocityFESpace.fe().nbFEDof(), nDimensions ),
493  M_elementDisplacement ( this->M_velocityFESpace.fe().nbFEDof(), nDimensions ),
494  M_elementVelocityRightHandSide ( this->M_velocityFESpace.fe().nbFEDof(), nDimensions ),
495  M_u_loc ( this->M_velocityFESpace.fe().nbFEDof(), nDimensions ),
496  M_reuseLinearPreconditioner ( true ),
497  M_mmFESpace ( &mmFESpace )
498 {
499 
500 }
501 
502 template<typename MeshType, typename SolverType>
503 OseenSolverShapeDerivative<MeshType, SolverType>::
504 ~OseenSolverShapeDerivative()
505 {
506 
507 }
508 
509 
510 // ===================================================
511 // Methods
512 // ===================================================
513 
514 template<typename MeshType, typename SolverType>
515 Real
516 OseenSolverShapeDerivative<MeshType, SolverType>::linearFlux ( const markerID_Type& flag )
517 {
518  return this->flux ( flag, M_linearSolution );
519 }
520 
521 template<typename MeshType, typename SolverType>
522 Real
523 OseenSolverShapeDerivative<MeshType, SolverType>::getLinearFlux ( const markerID_Type& flag )
524 {
525  if ( this->M_displayer->isLeader() )
526  {
527  std::cerr << "Warning: getLinearFlux is deprecated!" << std::endl
528  << " You should use linearFlux instead!" << std::endl;
529  }
530 
531  return linearFlux ( flag );
532 }
533 
534 template<typename MeshType, typename SolverType>
535 Real
536 OseenSolverShapeDerivative<MeshType, SolverType>::linearPressure ( const markerID_Type& flag )
537 {
538  return this->pressure ( flag, M_linearSolution );
539 }
540 
541 template<typename MeshType, typename SolverType>
542 Real
543 OseenSolverShapeDerivative<MeshType, SolverType>::getLinearPressure ( const markerID_Type& flag )
544 {
545  if ( this->M_displayer->isLeader() )
546  {
547  std::cerr << "Warning: getLinearPressure is deprecated!" << std::endl
548  << " You should use linearPressure instead!" << std::endl;
549  }
550 
551  return linearPressure ( flag );
552 }
553 
554 template<typename MeshType, typename SolverType>
555 Real
556 OseenSolverShapeDerivative<MeshType, SolverType>::linearLagrangeMultiplier ( const markerID_Type& flag, bcHandler_Type& bcHandler )
557 {
558  return lagrangeMultiplier ( flag, bcHandler, M_linearSolution );
559 }
560 
561 template<typename MeshType, typename SolverType>
562 Real
563 OseenSolverShapeDerivative<MeshType, SolverType>::getLinearLagrangeMultiplier ( const markerID_Type& flag, bcHandler_Type& bcHandler )
564 {
565  if ( this->M_displayer->isLeader() )
566  {
567  std::cerr << "Warning: getLinearLagrangeMultiplier is deprecated!" << std::endl
568  << " You should use linearLagrangeMultiplier instead!" << std::endl;
569  }
570 
571  return linearLagrangeMultiplier ( flag, bcHandler );
572 }
573 
574 template<typename MeshType, typename SolverType>
575 Real
576 OseenSolverShapeDerivative<MeshType, SolverType>::linearKineticNormalStress ( const markerID_Type& flag )
577 {
578  return this->linearKineticNormalStress ( flag, *this->M_solution, M_linearSolution );
579 }
580 
581 template<typename MeshType, typename SolverType>
582 Real
583 OseenSolverShapeDerivative<MeshType, SolverType>::linearKineticNormalStress ( const markerID_Type& flag, const vector_Type& solution, const vector_Type& linearSolution )
584 {
585  vector_Type velocityAndPressure ( solution, Repeated, Add );
586  vector_Type velocity ( this->M_velocityFESpace.map(), Repeated );
587  velocity.subset ( velocityAndPressure );
588 
589  vector_Type linearVelocityAndPressure ( linearSolution, Repeated );
590  vector_Type linearVelocity ( this->M_velocityFESpace.map(), Repeated );
591  linearVelocity.subset ( linearVelocityAndPressure );
592 
593  return this->M_postProcessing->kineticNormalStressDerivative ( velocity, linearVelocity, this->M_oseenData->density(), flag );
594 }
595 
596 template<typename MeshType, typename SolverType>
597 Real
598 OseenSolverShapeDerivative<MeshType, SolverType>::linearMeanNormalStress ( const markerID_Type& flag, bcHandler_Type& bcHandler )
599 {
600  return this->linearMeanNormalStress ( flag, bcHandler, M_linearSolution );
601 }
602 
603 template<typename MeshType, typename SolverType>
604 Real
605 OseenSolverShapeDerivative<MeshType, SolverType>::linearMeanNormalStress ( const markerID_Type& flag, bcHandler_Type& bcHandler, const vector_Type& linearSolution )
606 {
607  return this->meanNormalStress ( flag, bcHandler, linearSolution );
608 }
609 
610 template<typename MeshType, typename SolverType>
611 Real
612 OseenSolverShapeDerivative<MeshType, SolverType>::linearMeanTotalNormalStress ( const markerID_Type& flag, bcHandler_Type& bcHandler )
613 {
614  return this->meanNormalStress ( flag, bcHandler, M_linearSolution ) - linearKineticNormalStress ( flag, *this->M_solution, M_linearSolution );
615 }
616 
617 template<typename MeshType, typename SolverType>
618 Real
619 OseenSolverShapeDerivative<MeshType, SolverType>::linearMeanTotalNormalStress ( const markerID_Type& flag, bcHandler_Type& bcHandler, const vector_Type& solution, const vector_Type& linearSolution )
620 {
621  return this->meanNormalStress ( flag, bcHandler, linearSolution ) - linearKineticNormalStress ( flag, solution, linearSolution );
622 }
623 
624 template<typename MeshType, typename SolverType>
625 void OseenSolverShapeDerivative<MeshType, SolverType>::setUp ( const GetPot& dataFile )
626 {
627  // M_linearLinSolver.setDataFromGetPot( dataFile, "lin_fluid/solver" );
628  // M_linearLinSolver.setAztecooPreconditioner( dataFile, "lin_fluid/solver" );
629 
630  oseenSolver_Type::setUp ( dataFile );
631 
632  M_reuseLinearPreconditioner = dataFile ( "lin_fluid/prec/reuse", true );
633 
634  // std::string preconditionerType = dataFile( "lin_fluid/prec/prectype", "Ifpack" );
635 
636  // M_linearPreconditioner = prec_type( PRECFactory::instance().createObject( preconditionerType ) ); //linPrec is not used
637  // M_linearPreconditioner->setDataFromGetPot( dataFile, "lin_fluid/prec" );
638 
639 }
640 
641 
642 template<typename MeshType, typename SolverType>
643 void OseenSolverShapeDerivative<MeshType, SolverType>::solveLinearSystem ( bcHandler_Type& bcHandler )
644 {
645  this->M_Displayer.leaderPrint ( " LF- Finalizing the matrix and vectors ... " );
646 
647  LifeChrono chrono;
648  chrono.start();
649 
650  // matrix and vector assembling communication
651  this->M_matrixNoBC->globalAssemble();
652 
653  M_linearRightHandSideNoBC.globalAssemble();
654 
655  matrixPtr_Type matrixFull ( new matrix_Type ( this->M_localMap, this->M_matrixNoBC->meanNumEntries() ) );
656 
657  this->updateStabilization ( *matrixFull );
658  this->getFluidMatrix ( *matrixFull );
659 
660  vector_Type rightHandSideFull ( M_linearRightHandSideNoBC );
661 
662  chrono.stop();
663  this->M_Displayer.leaderPrintMax ( "done in " , chrono.diff() );
664 
665  // boundary conditions update
666  this->M_Displayer.leaderPrint ( " LF- Applying boundary conditions ... " );
667  chrono.start();
668 
669  this->applyBoundaryConditions ( *matrixFull, rightHandSideFull, bcHandler );
670 
671  chrono.stop();
672  this->M_Displayer.leaderPrintMax ( "done in ", chrono.diff() );
673 
674  // solving the system
675 
676  // using the same preconditioner as for the non linear problem (the matrix changes only in the
677  // boundary terms).
678  matrixFull->globalAssemble();
679  this->M_linearSolver->setMatrix ( *matrixFull );
680  this->M_linearSolver->setReusePreconditioner ( M_reuseLinearPreconditioner );
681  std::shared_ptr<MatrixEpetra<Real> > staticCast = std::static_pointer_cast<MatrixEpetra<Real> > (matrixFull);
682  this->M_linearSolver->solveSystem ( rightHandSideFull, M_linearSolution, staticCast );
683 
684  *this->M_residual = M_linearRightHandSideNoBC;
685  *this->M_residual -= *this->M_matrixNoBC * this->M_linearSolution;
686 
687  // if(S_verbose)
688  // {
689  // this->M_Displayer.leaderPrintMax( "NormInf Residual Lin = " , this->M_residual.NormInf());
690  // this->M_Displayer.leaderPrintMax( "NormInf Solution Lin = " , this->M_linearSolution.NormInf());
691  // }
692 } // solveLinearSystem
693 
694 template<typename MeshType, typename SolverType>
695 void
696 OseenSolverShapeDerivative<MeshType, SolverType>::updateLinearSystem ( const matrix_Type& /*matrixNoBC*/,
697  Real& /*alpha*/,
698  const vector_Type& un,
699  const vector_Type& uk,
700  const vector_Type& disp,
701  const vector_Type& w,
702  const vector_Type& dw,
703  const vector_Type& sourceVector )
704 {
705  this->M_Displayer.leaderPrint ( " LF- Updating the right hand side ... " );
706  LifeChrono chrono;
707  chrono.start();
708 
709  Int numVelocityComponent = nDimensions;
710  // sourceVector is usually zero
711  M_linearRightHandSideNoBC = sourceVector;
712 
713  if ( this->M_oseenData->useShapeDerivatives() )
714  {
715  // right hand side for the linearized ale system
716  // Loop on elements
717  vector_Type unRepeated ( un , Repeated );
718  vector_Type ukRepeated ( uk , Repeated );
719  vector_Type dispRepeated ( disp, Repeated );
720  vector_Type wRepeated ( w , Repeated );
721  vector_Type dwRepeated ( dw , Repeated );
722 
723  // std::cout << wRepeated.NormInf() << std::endl;
724  // std::cout << dwRepeated.NormInf() << std::endl;
725  // std::cout << dispRepeated.NormInf() << std::endl;
726 
727  vector_Type linearRightHandSideNoBC ( M_linearRightHandSideNoBC.map(), Repeated );
728 
729  for ( UInt i = 0; i < this->M_velocityFESpace.mesh()->numVolumes(); i++ )
730  {
731 
732  this->M_pressureFESpace.fe().update ( this->M_pressureFESpace.mesh()->volumeList ( i ) );
733  this->M_velocityFESpace.fe().updateFirstDerivQuadPt ( this->M_velocityFESpace.mesh()->volumeList ( i ) );
734 
735  // initialization of elementary vectors
736  M_elementVectorVelocity.zero();
737  M_elementVectorPressure.zero();
738 
739  for ( UInt iNode = 0 ; iNode < this->M_velocityFESpace.fe().nbFEDof() ; iNode++ )
740  {
741  UInt iLocal = this->M_velocityFESpace.fe().patternFirst ( iNode ); // iLocal = iNode
742 
743  for ( Int iComponent = 0; iComponent < numVelocityComponent; ++iComponent )
744  {
745  UInt iGlobal = this->M_velocityFESpace.dof().localToGlobalMap ( i, iLocal ) + iComponent * this->dimVelocity();
746 
747  // u^n - w^iNode local
748  M_elementConvectionVelocity.vec( ) [ iLocal + iComponent * this->M_velocityFESpace.fe().nbFEDof() ] = unRepeated (iGlobal)
749  - wRepeated ( iGlobal );
750  // w^iNode local
751  M_elementMeshVelocity.vec( ) [ iLocal + iComponent * this->M_velocityFESpace.fe().nbFEDof() ] = wRepeated ( iGlobal );
752  // u^iNode local
753  M_elementVelocity.vec( ) [ iLocal + iComponent * this->M_velocityFESpace.fe().nbFEDof() ] = ukRepeated ( iGlobal );
754  // d local
755  M_elementDisplacement.vec( ) [ iLocal + iComponent * this->M_velocityFESpace.fe().nbFEDof() ] = dispRepeated ( iGlobal );
756  // dw local
757  M_elementVelocityRightHandSide.vec( ) [ iLocal + iComponent * this->M_velocityFESpace.fe().nbFEDof() ] = dwRepeated ( iGlobal );
758  // un local
759  M_u_loc.vec() [ iLocal + iComponent * this->M_velocityFESpace.fe().nbFEDof() ] = unRepeated ( iGlobal );
760  }
761  }
762  /*
763  std::cout << M_elementConvectionVelocity.vec() << std::endl;
764  std::cout << M_elementMeshVelocity.vec() << std::endl;
765  std::cout << M_elementVelocity.vec() << std::endl;
766  std::cout << M_elementDisplacement.vec() << std::endl;
767  std::cout << M_elementVelocityRightHandSide.vec() << std::endl;
768  std::cout << M_u_loc.vec() << std::endl;
769  */
770  for ( UInt iNode = 0 ; iNode < this->M_pressureFESpace.fe().nbFEDof() ; iNode++ )
771  {
772  UInt iLocal = this->M_pressureFESpace.fe().patternFirst ( iNode ); // iLocal = iNode
773  UInt iGlobal = this->M_pressureFESpace.dof().localToGlobalMap ( i, iLocal ) + numVelocityComponent * this->dimVelocity();
774  M_elementPressure[ iLocal ] = ukRepeated[ iGlobal ]; // p^iNode local
775 
776  }
777 
778  // Elementary vectors
779 
780  //commented the code to print out the elementary data. Useful for debugging.
781 
782  // - \rho ( \grad( u^n-w^iNode ):[I\div d - (\grad d)^T] u^iNode + ( u^n-w^iNode )^T[I\div d - (\grad d)^T] (\grad u^iNode)^T , v )
783  AssemblyElemental::source_mass1 ( - this->M_oseenData->density(),
784  M_elementVelocity,
785  M_elementMeshVelocity,
786  M_elementConvectionVelocity,
787  M_elementDisplacement,
788  M_elementVectorVelocity,
789  this->M_velocityFESpace.fe() );
790  /*
791  std::cout << "source_mass1 -> norm_inf(M_elementVectorVelocity)" << std::endl;
792  M_elementVectorVelocity.showMe(std::cout);
793  */
794  // + \rho * ( \grad u^iNode dw, v )
795  AssemblyElemental::source_mass2 ( this->M_oseenData->density(),
796  M_elementVelocity,
797  M_elementVelocityRightHandSide,
798  M_elementVectorVelocity,
799  this->M_velocityFESpace.fe() );
800  /*
801  std::cout << "source_mass2 -> norm_inf(M_elementVectorVelocity)" << std::endl;
802  M_elementVectorVelocity.showMe(std::cout);
803  */
804  // - \rho/2 ( \nabla u^n:[2 * I\div d - (\grad d)^T] u^iNode , v )
805  AssemblyElemental::source_mass3 ( - 0.5 * this->M_oseenData->density(),
806  M_u_loc,
807  M_elementVelocity,
808  M_elementDisplacement,
809  M_elementVectorVelocity,
810  this->M_velocityFESpace.fe() );
811  /*
812  std::cout << "source_mass3 -> norm_inf(M_elementVectorVelocity)" << std::endl;
813  M_elementVectorVelocity.showMe(std::cout);
814  */
815  // - ( [-p^iNode I + 2*mu e(u^iNode)] [I\div d - (\grad d)^T] , \grad v )
816  AssemblyElemental::source_stress ( - 1.0,
817  this->M_oseenData->viscosity(),
818  M_elementVelocity,
819  M_elementPressure,
820  M_elementDisplacement,
821  M_elementVectorVelocity,
822  this->M_velocityFESpace.fe(),
823  this->M_pressureFESpace.fe() );
824  /*
825  std::cout << "source_stress -> norm_inf(M_elementVectorVelocity)" << std::endl;
826  M_elementVectorVelocity.showMe(std::cout);
827  */
828  // + \mu ( \grad u^iNode \grad d + [\grad d]^T[\grad u^iNode]^T : \grad v )
829  AssemblyElemental::source_stress2 ( this->M_oseenData->viscosity(),
830  M_elementVelocity,
831  M_elementDisplacement,
832  M_elementVectorVelocity,
833  this->M_velocityFESpace.fe() );
834  /*
835  std::cout << "source_stress2 -> norm_inf(M_elementVectorVelocity)" << std::endl;
836  M_elementVectorVelocity.showMe(std::cout);
837  */
838  // + ( (\grad u^iNode):[I\div d - (\grad d)^T] , q )
839  AssemblyElemental::source_press ( -1.0,
840  M_elementVelocity,
841  M_elementDisplacement,
842  M_elementVectorPressure,
843  this->M_velocityFESpace.fe(),
844  this->M_pressureFESpace.fe() );
845  /*
846  std::cout << "source_press -> norm_inf(M_elementVectorVelocity)" << std::endl;
847  M_elementVectorPressure.showMe(std::cout);
848  */
849  //
850  // Assembling
851  //
852  /*
853  std::cout << "debut ====================" << std::endl;
854  M_elementVectorPressure.showMe(std::cout);
855  M_elementVectorVelocity.showMe(std::cout);
856  std::cout << "fin ====================" << std::endl;
857  */
858  // assembling pressure right hand side
859  assembleVector ( linearRightHandSideNoBC,
860  M_elementVectorPressure,
861  this->M_pressureFESpace.fe(),
862  this->M_pressureFESpace.dof(),
863  0,
864  numVelocityComponent * this->dimVelocity() );
865 
866  // loop on velocity components
867  for ( Int iComponent = 0; iComponent < numVelocityComponent; iComponent++ )
868  {
869  // assembling velocity right hand side
870  assembleVector ( linearRightHandSideNoBC,
871  M_elementVectorVelocity,
872  this->M_velocityFESpace.fe(),
873  this->M_velocityFESpace.dof(),
874  iComponent,
875  iComponent * this->dimVelocity() );
876  }
877  }
878 
879  linearRightHandSideNoBC.globalAssemble();
880  M_linearRightHandSideNoBC += linearRightHandSideNoBC;
881  // M_linearRightHandSideNoBC *= -1.;
882  // if( S_verbose )
883  // this->M_Displayer.leaderPrint( "norm( M_linearRightHandSideNoBC) = " , M_linearRightHandSideNoBC.NormInf() );
884 
885  }
886 
887  chrono.stop();
888  this->M_Displayer.leaderPrintMax ( "done in ", chrono.diff() );
889 }
890 
891 
892 //#if UNDEF
893 template<typename MeshType, typename SolverType>
894 void
895 OseenSolverShapeDerivative<MeshType, SolverType>::
896 updateShapeDerivatives ( matrix_Type& matrix,
897  Real& alpha,
898  const vector_Type& un,
899  const vector_Type& uk,
900  //const vector_Type& disp,
901  const vector_Type& w,
902  UInt offset,
903  FESpace<mesh_Type, MapEpetra>& mmFESpace,
904  bool wImplicit,
905  bool convectiveTermDerivative )
906 {
907  LifeChrono chrono;
908 
909  UInt numVelocityComponent = nDimensions;
910 
911  // M_linearRightHandSideNoBC = sourceVector;//which is usually zero
912 
913  if ( this->M_oseenData->useShapeDerivatives() )
914  {
915  this->M_Displayer.leaderPrint ( " LF- Updating shape derivative blocks ... " );
916 
917  //
918  // RIGHT HAND SIDE FOR THE LINEARIZED ALE SYSTEM
919  //
920  chrono.start();
921 
922  // Loop on elements
923 
924  vector_Type unRepeated ( un , Repeated );
925  vector_Type ukRepeated ( uk , Repeated );
926  // vector_Type dispRepeated( disp, Repeated );
927  vector_Type wRepeated ( w , Repeated );
928  // vector_Type dwRepeated ( dw , Repeated );
929 
930  // std::cout << wRepeated.NormInf() << std::endl;
931  // std::cout << dwRepeated.NormInf() << std::endl;
932  // std::cout << dispRepeated.NormInf() << std::endl;
933 
934  // vector_Type rhsLinNoBC( M_linearRightHandSideNoBC.map(), Repeated);
935 
936  for ( UInt i = 0; i < this->M_velocityFESpace.mesh()->numVolumes(); i++ )
937  {
938 
939  this->M_pressureFESpace.fe().update ( this->M_pressureFESpace.mesh()->volumeList ( i ) );
940  this->M_velocityFESpace.fe().updateFirstDerivQuadPt ( this->M_velocityFESpace.mesh()->volumeList ( i ) );
941  this->M_pressureFESpace.fe().updateFirstDerivQuadPt ( this->M_velocityFESpace.mesh()->volumeList ( i ) );
942 
943  // just to provide the id number in the assem_mat_mixed
944  //this->M_pressureFESpace.fe().updateFirstDeriv( this->M_velocityFESpace.mesh()->volumeList( i ) );
945  //as updateFirstDer
946  //this->M_velocityFESpace.fe().updateFirstDeriv( this->M_velocityFESpace.mesh()->volumeList( i ) );
947  mmFESpace.fe().updateFirstDerivQuadPt ( mmFESpace.mesh()->volumeList ( i ) );
948 
949  // initialization of elementary vectors
950  std::shared_ptr<MatrixElemental> elementMatrixPressure ( new MatrixElemental ( this->M_pressureFESpace.fe().nbFEDof(),
951  1,
952  0,
953  mmFESpace.fe().nbFEDof(),
954  0,
955  nDimensions ) );
956  std::shared_ptr<MatrixElemental> elementMatrixVelocity ( new MatrixElemental ( this->M_velocityFESpace.fe().nbFEDof(),
957  nDimensions,
958  0,
959  this->M_velocityFESpace.fe().nbFEDof(),
960  0,
961  nDimensions ) );
962  std::shared_ptr<MatrixElemental> elementMatrixConvective;
963 
964  if ( convectiveTermDerivative )
965  {
966  elementMatrixConvective.reset ( new MatrixElemental ( this->M_velocityFESpace.fe().nbFEDof(),
967  nDimensions,
968  0,
969  this->M_velocityFESpace.fe().nbFEDof(),
970  0,
971  nDimensions ) );
972  elementMatrixConvective->zero();
973  }
974 
975  elementMatrixPressure->zero();
976  elementMatrixVelocity->zero();
977 
978  for ( UInt iNode = 0 ; iNode < this->M_velocityFESpace.fe().nbFEDof() ; iNode++ )
979  {
980  UInt iLocal = this->M_velocityFESpace.fe().patternFirst ( iNode ); // iLocal = iNode
981 
982  for ( UInt iComponent = 0; iComponent < numVelocityComponent; ++iComponent )
983  {
984  UInt iGlobal = this->M_velocityFESpace.dof().localToGlobalMap ( i, iLocal ) + iComponent * this->dimVelocity();
985 
986  // if(!wImplicit)
987  // u^n - w^iNode local
988  M_elementConvectionVelocity.vec() [ iLocal + iComponent * this->M_velocityFESpace.fe().nbFEDof() ] = unRepeated (iGlobal)
989  - wRepeated ( iGlobal );
990  // else
991  // u^n - w^iNode local
992  // M_elementConvectionVelocity.vec() [ iLocal + iComponent*this->M_velocityFESpace.fe().nbFEDof() ] = ukRepeated(iGlobal)
993  // - wRepeated(iGlobal);
994  // w^iNode local
995  M_elementMeshVelocity.vec( ) [ iLocal + iComponent * this->M_velocityFESpace.fe().nbFEDof() ] = wRepeated ( iGlobal );
996  // u^iNode local
997  M_elementVelocity.vec( ) [ iLocal + iComponent * this->M_velocityFESpace.fe().nbFEDof() ] = ukRepeated ( iGlobal );
998  // dw local
999  //M_elementDisplacement.vec( ) [ iLocal + iComponent*this->M_velocityFESpace.fe().nbFEDof() ] = dispRepeated( iGlobal );
1000  // dw local
1001  //M_elementVelocityRightHandSide.vec( ) [ iLocal + iComponent*this->M_velocityFESpace.fe().nbFEDof() ] = dwRepeated( iGlobal );
1002  // un local
1003  M_u_loc.vec() [ iLocal + iComponent * this->M_velocityFESpace.fe().nbFEDof() ] = unRepeated ( iGlobal );
1004  }
1005  }
1006  /*
1007  std::cout << M_elementConvectionVelocity.vec() << std::endl;
1008  std::cout << M_elementMeshVelocity.vec() << std::endl;
1009  std::cout << M_elementVelocity.vec() << std::endl;
1010  std::cout << M_elementDisplacement.vec() << std::endl;
1011  std::cout << M_elementVelocityRightHandSide.vec() << std::endl;
1012  std::cout << M_u_loc.vec() << std::endl;
1013  */
1014  for ( UInt iNode = 0 ; iNode < this->M_pressureFESpace.fe().nbFEDof() ; iNode++ )
1015  {
1016  // iLocal = iNode
1017  UInt iLocal = this->M_pressureFESpace.fe().patternFirst ( iNode );
1018  UInt iGlobal = this->M_pressureFESpace.dof().localToGlobalMap ( i, iLocal ) + numVelocityComponent * this->dimVelocity();
1019  // p^iNode local
1020  M_elementPressure[ iLocal ] = ukRepeated[ iGlobal ];
1021  }
1022 
1023  AssemblyElemental::shape_terms ( //M_elementDisplacement,
1024  this->M_oseenData->density(),
1025  this->M_oseenData->viscosity(),
1026  M_u_loc,
1027  M_elementVelocity,
1028  M_elementMeshVelocity,
1029  M_elementConvectionVelocity,
1030  M_elementPressure,
1031  *elementMatrixVelocity,
1032  mmFESpace.fe(),
1033  this->M_pressureFESpace.fe(),
1034  (ID) mmFESpace.fe().nbFEDof(),
1035  *elementMatrixPressure,
1036  0,
1037  wImplicit,
1038  alpha//,
1039  //elementMatrixConvective
1040  );
1041 
1042  //elementMatrixVelocity->showMe(std::cout);
1043 
1044  /*
1045  source_mass2( this->M_oseenData->density(),
1046  M_elementVelocity,
1047  *M_elementMatrixConvective,
1048  this->M_velocityFESpace.fe(),
1049  alpha );
1050  */
1051 
1052  AssemblyElemental::source_press ( 1.0,
1053  M_elementVelocity,
1054  *elementMatrixPressure,
1055  mmFESpace.fe(),
1056  this->M_pressureFESpace.fe(),
1057  (ID) mmFESpace.fe().nbFEDof() );
1058 
1059  //derivative of the convective term
1060  if ( convectiveTermDerivative )
1061  AssemblyElemental::mass_gradu ( this->M_oseenData->density(),
1062  M_elementVelocity,
1063  *elementMatrixConvective,
1064  this->M_velocityFESpace.fe() );
1065  /*
1066  std::cout << "source_press -> norm_inf( M_elementVectorVelocity )" << std::endl;
1067  M_elementVectorPressure.showMe( std::cout );
1068  */
1069  //
1070  // Assembling
1071  //
1072  /*
1073  std::cout << "debut ====================" << std::endl;
1074  M_elementVectorPressure.showMe( std::cout );
1075  M_elementVectorVelocity.showMe( std::cout );
1076  std::cout << "fin ====================" << std::endl;
1077  */
1078  UInt const velocityTotalDof ( this->M_velocityFESpace.dof().numTotalDof() );
1079  UInt const meshTotalDof ( mmFESpace.dof().numTotalDof() );
1080  for ( UInt iComponent = 0; iComponent < numVelocityComponent; ++iComponent )
1081  {
1082  for ( UInt jComponent = 0; jComponent < numVelocityComponent; ++jComponent )
1083  {
1084  assembleMatrix ( matrix,
1085  *elementMatrixVelocity,
1086  this->M_velocityFESpace.fe(),
1087  mmFESpace.fe(),
1088  this->M_velocityFESpace.dof(),
1089  mmFESpace.dof(),
1090  iComponent,
1091  jComponent,
1092  iComponent * velocityTotalDof,
1093  offset + jComponent * meshTotalDof );
1094 
1095  //assembling the derivative of the convective term
1096  if ( convectiveTermDerivative )
1097  assembleMatrix ( matrix,
1098  *elementMatrixConvective,
1099  this->M_velocityFESpace.fe(),
1100  this->M_velocityFESpace.fe(),
1101  this->M_velocityFESpace.dof(),
1102  this->M_velocityFESpace.dof(),
1103  iComponent,
1104  jComponent,
1105  iComponent * velocityTotalDof,
1106  jComponent * velocityTotalDof );
1107  }
1108 
1109  assembleMatrix ( matrix,
1110  *elementMatrixPressure,
1111  this->M_pressureFESpace.fe(),
1112  mmFESpace.fe(),
1113  this->M_pressureFESpace.dof(),
1114  mmFESpace.dof(),
1115  (UInt) 0,
1116  iComponent,
1117  (UInt) numVelocityComponent * velocityTotalDof,
1118  offset + iComponent * meshTotalDof );
1119  }
1120  }
1121  }
1122 
1123  chrono.stop();
1124  this->M_Displayer.leaderPrintMax ("done in ", chrono.diff() );
1125 }
1126 
1127 } // namespace LifeV
1128 
1129 #endif // OSEENSOLVERSHAPEDERIVATIVE_H
LifeV::LifeChrono::start
void start()
Start the timer.
Definition: LifeChrono.hpp:93
LifeV::SolverAztecOO
SolverAztecOO - Class to wrap linear solver.
Definition: SolverAztecOO.hpp:78
LifeV::Int
int32_type Int
Generic integer data.
Definition: LifeV.hpp:188
LifeV::markerID_Type
ID markerID_Type
markerID_Type is the type used to store the geometric entity marker IDs
Definition: Marker.hpp:81
ETCurrentFE::updateInverseJacobian
void updateInverseJacobian(const UInt &iQuadPt)
Definition: ETCurrentFE.cpp:405
LifeV::MapEpetra::importer
Epetra_Import const & importer()
Getter for the Epetra_Import.
Definition: MapEpetra.cpp:394
LifeV::LifeChrono::diff
Real diff()
Compute the difference in time between start and stop.
Definition: LifeChrono.hpp:111
LifeV::Real
double Real
Generic real data.
Definition: LifeV.hpp:175
LifeV::nDimensions
const UInt nDimensions(NDIM)
LifeV::LifeChrono::stop
void stop()
Stop the timer.
Definition: LifeChrono.hpp:100
LIFEV_DEPRECATED
#define LIFEV_DEPRECATED(func)
Definition: LifeV.hpp:117
GetPot::operator()
bool operator()(const char *VarName, bool Default) const
Definition: GetPot.hpp:2008
LifeV::UInt
uint32_type UInt
generic unsigned integer (used mainly for addressing)
Definition: LifeV.hpp:191