OOFEM  2.4
OOFEM.org - Object Oriented Finite Element Solver
nlinearstatic.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 "../sm/EngineeringModels/nlinearstatic.h"
36 #include "../sm/Elements/structuralelement.h"
37 #include "nummet.h"
38 #include "timestep.h"
39 #include "metastep.h"
40 #include "error.h"
41 #include "verbose.h"
42 #include "sparsenonlinsystemnm.h"
43 #include "nrsolver.h"
44 #include "calmls.h"
45 #include "outputmanager.h"
46 #include "datastream.h"
47 #include "classfactory.h"
48 #include "timer.h"
49 #include "contextioerr.h"
50 #include "sparsemtrx.h"
51 #include "errorestimator.h"
52 #include "mathfem.h"
53 #include "dofmanager.h"
54 #include "dof.h"
55 #include "unknownnumberingscheme.h"
56 
57 #ifdef __PARALLEL_MODE
58  #include "problemcomm.h"
59  #include "communicator.h"
60  #include "loadbalancer.h"
61 #endif
62 
63 namespace oofem {
64 REGISTER_EngngModel(NonLinearStatic);
65 
67  totalDisplacement(), incrementOfDisplacement(), internalForces(), initialLoadVector(), incrementalLoadVector(),
68  initialLoadVectorOfPrescribed(), incrementalLoadVectorOfPrescribed()
69 {
70  //
71  // constructor
72  //
73  prevStepLength = 0.;
74  currentStepLength = 0.;
78  stiffMode = nls_tangentStiffness; // default
80  initFlag = loadInitFlag = 1;
83  nMethod = NULL;
85 }
86 
87 
89 {
90  delete nMethod;
91 }
92 
93 
95 {
96  IRResultType result; // Required by IR_GIVE_FIELD macro
97 
98  if ( mStep == NULL ) {
99  OOFEM_ERROR("undefined meta step");
100  }
101 
102  int _val = 0;
104  IR_GIVE_OPTIONAL_FIELD( ( mStep->giveAttributesRecord() ), _val, "controllmode" );
106 
107  if ( mode == nls_indirectControl ) {
108  if ( nMethod ) {
109  if ( dynamic_cast< CylindricalALM * >(nMethod) ) {
110  return nMethod;
111  } else {
112  delete nMethod;
113  }
114  }
115 
116  this->nMethod = new CylindricalALM(this->giveDomain(1), this);
117  } else if ( mode == nls_directControl ) {
118  if ( nMethod ) {
119  if ( dynamic_cast< NRSolver * >(nMethod) ) {
120  return nMethod;
121  } else {
122  delete nMethod;
123  }
124  }
125 
126  this->nMethod = new NRSolver(this->giveDomain(1), this);
127  } else {
128  OOFEM_ERROR("unsupported controlMode");
129  }
130 
131  return this->nMethod;
132 }
133 
134 
135 void
137 {
138  IRResultType result; // Required by IR_GIVE_FIELD macro
139 
140  MetaStep *mStep1 = this->giveMetaStep( mStep->giveNumber() ); //this line ensures correct input file in staggered problem
141  InputRecord *ir = mStep1->giveAttributesRecord();
142 
144 
145  /*
146  * if ((mstep->giveFirstStepNumber() == tStep->giveNumber()) && hasString(initString, "fixload")) {
147  * double factor;
148  *
149  * printf ("NonLinearStatic: fixed load level");
150  * if (initialLoadVector.isEmpty()) initialLoadVector.resize(loadVector.giveSize());
151  * if ((controlMode == nls_directControl) || (controlMode == nls_directControl2)) factor = 1.0;
152  * else factor = loadLevel;
153  * loadVector.times (factor);
154  * initialLoadVector.add(loadVector);
155  * loadVector.zero();
156  * this->loadInitFlag = 1;
157  * this->loadLevel = 0.0;
158  * }
159  */
160  int _val = nls_indirectControl;
162  IR_GIVE_OPTIONAL_FIELD(ir, _val, "controllmode");
163  this->controlMode = ( NonLinearStatic_controlType ) _val;
164 
165  _val = IG_None;
167  this->initialGuessType = ( InitialGuess ) _val;
168 
169  deltaT = 1.0;
171  if ( deltaT < 0. ) {
172  OOFEM_ERROR("deltaT < 0");
173  }
174 
175  _val = nls_tangentStiffness;
177  this->stiffMode = ( NonLinearStatic_stiffnessMode ) _val;
178 
182 
184 
185  // called just to mark field as recognized, used later
186  ir->hasField(_IFT_NonLinearStatic_donotfixload);
187 }
188 
189 
192 {
193  IRResultType result; // Required by IR_GIVE_FIELD macro
194 
195  result = LinearStatic :: initializeFrom(ir);
196  if ( result != IRRT_OK ) {
197  return result;
198  }
199 
202 
206  }
207 
208 #ifdef __PARALLEL_MODE
209  if ( isParallel() ) {
210  //commBuff = new CommunicatorBuff (this->giveNumberOfProcesses(), CBT_dynamic);
212  communicator = new NodeCommunicator(this, commBuff, this->giveRank(),
213  this->giveNumberOfProcesses());
214 
216  nonlocalExt = 1;
218  this->giveNumberOfProcesses());
219  }
220  }
221 #endif
222 
223  return IRRT_OK;
224 }
225 
226 
228 // returns unknown quantity like displacement, velocity of equation eq
229 // This function translates this request to numerical method language
230 {
231  int eq = dof->__giveEquationNumber();
232 #ifdef DEBUG
233  if ( eq == 0 ) {
234  OOFEM_ERROR("invalid equation number");
235  }
236 #endif
237 
238 
239  if ( tStep != this->giveCurrentStep() ) {
240  OOFEM_ERROR("unknown time step encountered");
241  return 0.;
242  }
243 
244  switch ( mode ) {
245  case VM_Incremental:
246  // return incrementOfDisplacement -> at(eq);
247  // return nMethod-> giveUnknownComponent(IncrementOfSolution, eq);
249  return incrementOfDisplacement.at(eq);
250  } else {
251  return 0.;
252  }
253 
254  case VM_Total:
255  if ( totalDisplacement.isNotEmpty() ) {
256  return totalDisplacement.at(eq);
257  } else {
258  return 0.;
259  }
260 
261  case VM_Velocity:
263  return incrementOfDisplacement.at(eq) / tStep->giveTimeIncrement();
264  } else {
265  return 0.;
266  }
267 
268  default:
269  OOFEM_ERROR("Unknown is of undefined ValueModeType for this problem");
270  }
271 
272  return 0.0;
273 }
274 
276 {
277  if ( master && (!force)) {
279  } else {
280  if ( !stepWhenIcApply ) {
281  int inin = giveNumberOfTimeStepWhenIcApply();
282  // int nFirst = giveNumberOfFirstStep();
283  stepWhenIcApply.reset(new TimeStep(inin, this, 0, -deltaT, deltaT, 0));
284  }
285 
286  return stepWhenIcApply.get();
287  }
288 }
289 
290 
292 {
293  int istep = giveNumberOfFirstStep();
294  int mStepNum = 1;
295  double totalTime = 0.0;
296  StateCounterType counter = 1;
297  double deltaTtmp = deltaT;
298 
299  //do not increase deltaT on microproblem
300  if ( pScale == microScale ) {
301  deltaTtmp = 0.;
302  }
303 
304  if ( currentStep ) {
305  totalTime = currentStep->giveTargetTime() + deltaTtmp;
306  istep = currentStep->giveNumber() + 1;
307  counter = currentStep->giveSolutionStateCounter() + 1;
308  mStepNum = currentStep->giveMetaStepNumber();
309 
310  if ( !this->giveMetaStep(mStepNum)->isStepValid(istep) ) {
311  mStepNum++;
312  if ( mStepNum > nMetaSteps ) {
313  OOFEM_ERROR("no next step available, mStepNum=%d > nMetaSteps=%d", mStepNum, nMetaSteps);
314  }
315  }
316  } else {
317  // first step -> generate initial step
319  currentStep.reset(new TimeStep(*newStep));
320  }
321 
322  previousStep = std :: move(currentStep);
323  currentStep.reset( new TimeStep(istep, this, mStepNum, totalTime, deltaTtmp, counter) );
324  // dt variable are set eq to 0 for statics - has no meaning
325  // *Wrong* It has meaning for viscoelastic materials.
326 
327  return currentStep.get();
328 }
329 
330 
332 {
333  if ( this->isParallel() ) {
334  #ifdef __VERBOSE_PARALLEL
335  // force equation numbering before setting up comm maps
337  OOFEM_LOG_INFO("[process rank %d] neq is %d\n", this->giveRank(), neq);
338  #endif
339 
340  // set up communication patterns
341  this->initializeCommMaps();
342  // init remote dofman list
343  // this->initRemoteDofManList ();
344  }
346 }
347 
348 
349 void
351 {
352  proceedStep(1, tStep);
353 }
354 
355 
356 void
358 {
359  this->doStepOutput(tStep);
360  this->updateLoadVectors(tStep);
361  this->saveStepContext(tStep, CM_State | CM_Definition);
362 }
363 
364 
365 void
367 {
368  MetaStep *mstep = this->giveMetaStep( tStep->giveMetaStepNumber() );
369  bool isLastMetaStep = ( tStep->giveNumber() == mstep->giveLastStepNumber() );
370 
371  if ( controlMode == nls_indirectControl ) {
372  //if ((tStep->giveNumber() == mstep->giveLastStepNumber()) && ir->hasField("fixload")) {
373  if ( isLastMetaStep ) {
375  OOFEM_LOG_INFO("Fixed load level\n");
376 
377  //update initialLoadVector
379 
381 
384 
385  this->loadInitFlag = 1;
386  }
387 
388  //if (!mstep->giveAttributesRecord()->hasField("keepll")) this->loadLevelInitFlag = 1;
389  }
390  } else { // direct control
391  //update initialLoadVector after each step of direct control
392  //(here the loading is not proportional)
393  OOFEM_LOG_DEBUG("Fixed load level\n");
394 
396 
398 
401 
402  this->loadInitFlag = 1;
403  }
404 
405 
406  // if (isLastMetaStep) {
407  if ( isLastMetaStep && !mstep->giveAttributesRecord()->hasField(_IFT_NonLinearStatic_donotfixload) ) {
408 #ifdef VERBOSE
409  OOFEM_LOG_INFO("Reseting load level\n");
410 #endif
413  } else {
414  cumulatedLoadLevel = 0.0;
415  }
416 
417  this->loadLevel = 0.0;
418  }
419 }
420 
421 
422 void
424 {
425  //
426  // creates system of governing eq's and solves them at given time step
427  //
428  // first assemble problem at current time step
429  //
430 
431  if ( initFlag ) {
432  //
433  // first step create space for stiffness Matrix
434  //
435  if ( !stiffnessMatrix ) {
437  }
438 
439  if ( !stiffnessMatrix ) {
440  OOFEM_ERROR("sparse matrix creation failed");
441  }
442 
443  if ( nonlocalStiffnessFlag ) {
444  if ( !stiffnessMatrix->isAsymmetric() ) {
445  OOFEM_ERROR("stiffnessMatrix does not support asymmetric storage");
446  }
447  }
448 
449  stiffnessMatrix->buildInternalStructure( this, di, EModelDefaultEquationNumbering() );
450  }
451 
452 #if 0
453  if ( ( mstep->giveFirstStepNumber() == tStep->giveNumber() ) ) {
454  #ifdef VERBOSE
455  OOFEM_LOG_INFO("Resetting load level\n");
456  #endif
459  } else {
460  cumulatedLoadLevel = 0.0;
461  }
462  this->loadLevel = 0.0;
463  }
464 #endif
465 
467 #ifdef VERBOSE
468  OOFEM_LOG_DEBUG("Assembling reference load\n");
469 #endif
470  //
471  // assemble the incremental reference load vector
472  //
474  refLoadInputMode, this->giveDomain(di), tStep);
475 
476  loadInitFlag = 0;
477  }
478 
479  if ( tStep->giveNumber() == 1 ) {
485  }
486 
487  //
488  // -> BEGINNING OF LOAD (OR DISPLACEMENT) STEP <-
489  //
490 
491  //
492  // set-up numerical model
493  //
494  this->giveNumericalMethod( this->giveMetaStep( tStep->giveMetaStepNumber() ) );
495  //
496  // call numerical model to solve arise problem
497  //
498 #ifdef VERBOSE
499  OOFEM_LOG_RELEVANT( "\n\nSolving [step number %5d.%d, time = %e]\n\n", tStep->giveNumber(), tStep->giveVersion(), tStep->giveIntrinsicTime() );
500 #endif
501 
502  if ( this->initialGuessType == IG_Tangent ) {
503 #ifdef VERBOSE
504  OOFEM_LOG_RELEVANT("Computing initial guess\n");
505 #endif
506  FloatArray extrapolatedForces;
507  this->assemblePrescribedExtrapolatedForces( extrapolatedForces, tStep, TangentStiffnessMatrix, this->giveDomain(di) );
508  extrapolatedForces.negated();
509  this->updateComponent( tStep, NonLinearLhs, this->giveDomain(di) );
511  OOFEM_LOG_RELEVANT("solving for increment\n");
512  linSolver->solve(*stiffnessMatrix, extrapolatedForces, incrementOfDisplacement);
513  OOFEM_LOG_RELEVANT("initial guess found\n");
515  } else if ( this->initialGuessType != IG_None ) {
516  OOFEM_ERROR("Initial guess type: %d not supported", initialGuessType);
517  } else {
519  }
520 
521  //totalDisplacement.printYourself();
522  if ( initialLoadVector.isNotEmpty() ) {
526  } else {
530  }
532  //this->updateComponent(tStep, NonLinearLhs, this->giveDomain(di));
534  OOFEM_LOG_RELEVANT("Equilibrium reached at load level = %f in %d iterations\n", cumulatedLoadLevel + loadLevel, currentIterations);
536 }
537 
538 void
540 //
541 // updates some component, which is used by numerical method
542 // to newly reached state. used mainly by numerical method
543 // when new tangent stiffness is needed during finding
544 // of new equilibrium stage.
545 //
546 {
547  switch ( cmpn ) {
548  case NonLinearLhs:
549  if ( stiffMode == nls_tangentStiffness ) {
550  stiffnessMatrix->zero(); // zero stiffness matrix
551 #ifdef VERBOSE
552  OOFEM_LOG_DEBUG("Assembling tangent stiffness matrix\n");
553 #endif
554  this->assemble(*stiffnessMatrix, tStep, TangentAssembler(TangentStiffness),
557 #ifdef VERBOSE
558  OOFEM_LOG_DEBUG("Assembling secant stiffness matrix\n");
559 #endif
560  stiffnessMatrix->zero(); // zero stiffness matrix
561  this->assemble(*stiffnessMatrix, tStep, TangentAssembler(SecantStiffness),
563  initFlag = 0;
564  } else if ( ( stiffMode == nls_elasticStiffness ) && ( initFlag ||
565  ( this->giveMetaStep( tStep->giveMetaStepNumber() )->giveFirstStepNumber() == tStep->giveNumber() ) || (updateElasticStiffnessFlag) ) ) {
566 #ifdef VERBOSE
567  OOFEM_LOG_DEBUG("Assembling elastic stiffness matrix\n");
568 #endif
569  stiffnessMatrix->zero(); // zero stiffness matrix
570  this->assemble(*stiffnessMatrix, tStep, TangentAssembler(ElasticStiffness),
572  initFlag = 0;
573  } else {
574  // currently no action , this method is mainly intended to
575  // assemble new tangent stiffness after each iteration
576  // when secantStiffMode is on, we use the same stiffness
577  // during iteration process
578  }
579 
580  break;
581  case InternalRhs:
582 #ifdef VERBOSE
583  OOFEM_LOG_DEBUG("Updating internal forces\n");
584 #endif
585  // update internalForces and internalForcesEBENorm concurrently
586  this->giveInternalForces(internalForces, true, d->giveNumber(), tStep);
587  break;
588 
589  default:
590  OOFEM_ERROR("Unknown Type of component.");
591  }
592 }
593 
594 
595 void
597 {
598  if ( !this->giveDomain(1)->giveOutputManager()->testTimeStepOutput(tStep) ) {
599  return; // do not print even Solution step header
600  }
601 
602  fprintf(file, "\n\nOutput for time %.3e, solution step number %d\n", tStep->giveTargetTime(), tStep->giveNumber() );
603  fprintf(file, "Reached load level : %20.6f in %d iterations\n\n",
605 
606  nMethod->printState(file);
607 
608  this->giveDomain(1)->giveOutputManager()->doDofManOutput(file, tStep);
609  this->giveDomain(1)->giveOutputManager()->doElementOutput(file, tStep);
610  this->printReactionForces(tStep, 1, file);
611 }
612 
613 
616 {
617  contextIOResultType iores;
618 
619  if ( ( iores = EngngModel :: saveContext(stream, mode) ) != CIO_OK ) {
620  THROW_CIOERR(iores);
621  }
622 
623  if ( ( iores = totalDisplacement.storeYourself(stream) ) != CIO_OK ) {
624  THROW_CIOERR(iores);
625  }
626 
627  if ( ( iores = incrementOfDisplacement.storeYourself(stream) ) != CIO_OK ) {
628  THROW_CIOERR(iores);
629  }
630 
631  int _cm = controlMode;
632  if ( !stream.write(_cm) ) {
634  }
635 
636  if ( !stream.write(loadLevel) ) {
638  }
639 
640  if ( !stream.write(cumulatedLoadLevel) ) {
642  }
643 
644  if ( ( iores = initialLoadVector.storeYourself(stream) ) != CIO_OK ) {
645  THROW_CIOERR(iores);
646  }
647 
648  if ( ( iores = initialLoadVectorOfPrescribed.storeYourself(stream) ) != CIO_OK ) {
649  THROW_CIOERR(iores);
650  }
651 
652  return CIO_OK;
653 }
654 
655 
658 {
659  contextIOResultType iores;
660 
661  if ( ( iores = EngngModel :: restoreContext(stream, mode) ) != CIO_OK ) {
662  THROW_CIOERR(iores);
663  }
664 
665  //if ((iores = this->giveNumericalMethod(giveCurrentStep())->restoreContext (stream)) !=CIO_OK) THROW_CIOERR(iores);
666 
667  if ( ( iores = totalDisplacement.restoreYourself(stream) ) != CIO_OK ) {
668  THROW_CIOERR(iores);
669  }
670 
671  if ( ( iores = incrementOfDisplacement.restoreYourself(stream) ) != CIO_OK ) {
672  THROW_CIOERR(iores);
673  }
674 
675  int _cm;
676  if ( !stream.read(_cm) ) {
678  }
679 
681  if ( !stream.read(loadLevel) ) {
683  }
684 
685  if ( !stream.read(cumulatedLoadLevel) ) {
687  }
688 
689 
690  // store InitialLoadVector
691  if ( ( iores = initialLoadVector.restoreYourself(stream) ) != CIO_OK ) {
692  THROW_CIOERR(iores);
693  }
694 
695  if ( ( iores = initialLoadVectorOfPrescribed.restoreYourself(stream) ) != CIO_OK ) {
696  THROW_CIOERR(iores);
697  }
698 
699  return CIO_OK;
700 }
701 
702 
703 void
705 {
707 
708  this->giveNumericalMethod( this->giveCurrentMetaStep() )->setDomain( this->giveDomain(1) );
709 #ifdef __PARALLEL_MODE
710  if ( this->giveLoadBalancer() ) {
711  this->giveLoadBalancer()->setDomain( this->giveDomain(1) );
712  }
713 
714 #endif
715 }
716 
717 
718 void
720  const UnknownNumberingScheme &s, Domain *domain)
721 {
722 #ifdef TIME_REPORT
723  Timer timer;
724  timer.startTimer();
725 #endif
726 
727  LinearStatic :: assemble(answer, tStep, ma, s, domain);
728 
729  if ( ( nonlocalStiffnessFlag ) && dynamic_cast< const TangentAssembler* >(&ma) ) {
730  // add nonlocal contribution
731  for ( auto &elem : domain->giveElements() ) {
732  static_cast< StructuralElement * >( elem.get() )->addNonlocalStiffnessContributions(answer, s, tStep);
733  }
734 
735  // print storage statistics
736  answer.printStatistics();
737  }
738 
739 #ifdef TIME_REPORT
740  timer.stopTimer();
741  OOFEM_LOG_DEBUG( "NonLinearStatic: User time consumed by assembly: %.2fs\n", timer.getUtime() );
742 #endif
743 }
744 
745 
746 #ifdef __OOFEG
747 void
749 {
750  Domain *domain = this->giveDomain(1);
751  CharType ctype;
752 
753  if ( type != 1 ) {
754  return;
755  }
756 
757  if ( stiffMode == nls_tangentStiffness ) {
758  ctype = TangentStiffnessMatrix;
759  } else if ( stiffMode == nls_secantStiffness ) {
760  ctype = SecantStiffnessMatrix;
761  } else {
762  ctype = SecantStiffnessMatrix;
763  }
764 
765  for ( auto &elem : domain->giveElements() ) {
766  elem->showSparseMtrxStructure(ctype, gc, tStep);
767  }
768 
769  for ( auto &elem : domain->giveElements() ) {
770  elem->showExtendedSparseMtrxStructure(ctype, gc, tStep);
771  }
772 }
773 #endif
774 
775 
776 void
778 {
779  if ( ( di == 1 ) && ( tStep == this->giveCurrentStep() ) ) {
780  reactions = initialLoadVectorOfPrescribed;
782  } else {
783  OOFEM_ERROR("unable to respond due to invalid solution step or domain");
784  }
785 }
786 
787 
788 void
790  FloatArray &_incrementalLoadVectorOfPrescribed,
792  Domain *sourceDomain, TimeStep *tStep)
793 {
794  _incrementalLoadVector.resize( sourceDomain->giveEngngModel()->giveNumberOfDomainEquations( sourceDomain->giveNumber(), EModelDefaultEquationNumbering() ) );
795  _incrementalLoadVector.zero();
796  _incrementalLoadVectorOfPrescribed.resize( sourceDomain->giveEngngModel()->giveNumberOfDomainEquations( sourceDomain->giveNumber(), EModelDefaultPrescribedEquationNumbering() ) );
797  _incrementalLoadVectorOfPrescribed.zero();
798 
799  if ( _refMode == SparseNonLinearSystemNM :: rlm_incremental ) {
801  this->assembleVector(_incrementalLoadVector, tStep, ExternalForceAssembler(),
802  VM_Incremental, EModelDefaultEquationNumbering(), sourceDomain);
803 
804  this->assembleVector(_incrementalLoadVectorOfPrescribed, tStep, ExternalForceAssembler(),
805  VM_Incremental, EModelDefaultPrescribedEquationNumbering(), sourceDomain);
806  } else {
807  this->assembleVector(_incrementalLoadVector, tStep, ExternalForceAssembler(),
808  VM_Total, EModelDefaultEquationNumbering(), sourceDomain);
809 
810  this->assembleVector(_incrementalLoadVectorOfPrescribed, tStep, ExternalForceAssembler(),
811  VM_Total, EModelDefaultPrescribedEquationNumbering(), sourceDomain);
812  }
813 
815 }
816 
817 
818 int
819 NonLinearStatic :: estimateMaxPackSize(IntArray &commMap, DataStream &buff, int packUnpackType)
820 {
821  int count = 0, pcount = 0;
822  Domain *domain = this->giveDomain(1);
823 
824  if ( packUnpackType == 0 ) {
825  for ( int map: commMap ) {
826  DofManager *dman = domain->giveDofManager( map );
827  for ( Dof *dof: *dman ) {
828  if ( dof->isPrimaryDof() && ( dof->__giveEquationNumber() ) ) {
829  count++;
830  } else {
831  pcount++;
832  }
833  }
834  }
835 
836  //printf ("\nestimated count is %d\n",count);
837  return ( buff.givePackSizeOfDouble(1) * max(count, pcount) );
838  } else if ( packUnpackType == 1 ) {
839  for ( int map: commMap ) {
840  count += domain->giveElement( map )->estimatePackSize(buff);
841  }
842 
843  return count;
844  }
845 
846  return 0;
847 }
848 
849 
850 #ifdef __PARALLEL_MODE
851 LoadBalancer *
853 {
854  if ( lb ) {
855  return lb;
856  }
857 
858  if ( loadBalancingFlag ) {
860  lb = classFactory.createLoadBalancer( "parmetis", this->giveDomain(1) );
861  return lb;
862  } else {
863  return NULL;
864  }
865 }
866 
867 
870 {
871  if ( lbm ) {
872  return lbm;
873  }
874 
875  if ( loadBalancingFlag ) {
876  lbm = classFactory.createLoadBalancerMonitor( "wallclock", this);
877  return lbm;
878  } else {
879  return NULL;
880  }
881 }
882 #endif
883 
884 void
886 {
887  Domain *domain = this->giveDomain(1);
888  int ndofman = domain->giveNumberOfDofManagers();
889  bool initialLoadVectorEmpty = initialLoadVector.isEmpty();
890  bool initialLoadVectorOfPrescribedEmpty = initialLoadVectorOfPrescribed.isEmpty();
891 
892  for ( int idofman = 1; idofman <= ndofman; idofman++ ) {
893  DofManager *_dm = domain->giveDofManager(idofman);
894  for ( Dof *_dof: *_dm ) {
895  if ( _dof->isPrimaryDof() ) {
896  int _eq;
897  if ( ( _eq = _dof->__giveEquationNumber() ) ) {
898  // pack values in solution vectors
899  _dof->updateUnknownsDictionary( tStep, VM_Total, totalDisplacement.at(_eq) );
900  if ( initialLoadVectorEmpty ) {
901  _dof->updateUnknownsDictionary(tStep, VM_RhsInitial, 0.0);
902  } else {
903  _dof->updateUnknownsDictionary( tStep, VM_RhsInitial, initialLoadVector.at(_eq) );
904  }
905 
906  _dof->updateUnknownsDictionary( tStep, VM_RhsIncremental, incrementalLoadVector.at(_eq) );
907  } else if ( ( _eq = _dof->__givePrescribedEquationNumber() ) ) {
908  // pack values in prescribed solution vectors
909  if ( initialLoadVectorOfPrescribedEmpty ) {
910  _dof->updateUnknownsDictionary(tStep, VM_RhsInitial, 0.0);
911  } else {
912  _dof->updateUnknownsDictionary( tStep, VM_RhsInitial, initialLoadVectorOfPrescribed.at(_eq) );
913  }
914 
915  _dof->updateUnknownsDictionary( tStep, VM_RhsIncremental, incrementalLoadVectorOfPrescribed.at(_eq) );
916  }
917  } // end primary dof
918  } // end dof loop
919  } // end dofman loop
920 }
921 
922 
923 void
925 {
926  Domain *domain = this->giveDomain(1);
927  int ndofman = domain->giveNumberOfDofManagers();
928  //int myrank = this->giveRank();
929 
930  // resize target arrays
938 
939  for ( int idofman = 1; idofman <= ndofman; idofman++ ) {
940  DofManager *_dm = domain->giveDofManager(idofman);
941  for ( Dof *_dof: *_dm ) {
942  if ( _dof->isPrimaryDof() ) {
943  int _eq;
944  if ( ( _eq = _dof->__giveEquationNumber() ) ) {
945  // pack values in solution vectors
946  totalDisplacement.at(_eq) = _dof->giveUnknownsDictionaryValue( tStep, VM_Total );
947  initialLoadVector.at(_eq) = _dof->giveUnknownsDictionaryValue( tStep, VM_RhsInitial );
948  incrementalLoadVector.at(_eq) = _dof->giveUnknownsDictionaryValue( tStep, VM_RhsIncremental );
949 
950  #if 0
951  // debug print
952  if ( _dm->giveParallelMode() == DofManager_shared ) {
953  fprintf(stderr, "[%d] Shared: %d(%d) -> %d\n", myrank, idofman, idof, _eq);
954  } else {
955  fprintf(stderr, "[%d] Local : %d(%d) -> %d\n", myrank, idofman, idof, _eq);
956  }
957 
958  #endif
959  } else if ( ( _eq = _dof->__givePrescribedEquationNumber() ) ) {
960  // pack values in prescribed solution vectors
961  initialLoadVectorOfPrescribed.at(_eq) = _dof->giveUnknownsDictionaryValue( tStep, VM_RhsInitial );
962  incrementalLoadVectorOfPrescribed.at(_eq) = _dof->giveUnknownsDictionaryValue( tStep, VM_RhsIncremental );
963 
964  #if 0
965  // debug print
966  fprintf(stderr, "[%d] %d(%d) -> %d\n", myrank, idofman, idof, -_eq);
967  #endif
968  }
969  } // end primary dof
970  } // end dof loop
971  } // end dofman loop
972 
973  this->initializeCommMaps(true);
975  // reinitialize error estimator (if any)
976  if ( this->giveDomainErrorEstimator(1) ) {
978  }
979 
980  initFlag = true;
981 }
982 
983 } // end namespace oofem
EngngModelTimer timer
E-model timer.
Definition: engngm.h:267
virtual void updateAttributes(MetaStep *mStep)
Update receiver attributes according to step metaStep attributes.
Definition: engngm.C:589
virtual SparseLinearSystemNM * giveLinearSolver()
Constructs (if necessary) and returns a linear solver.
The representation of EngngModel default unknown numbering.
Solves an approximated tangent problem from the last iteration. Useful for changing Dirichlet boundar...
Definition: engngm.h:199
virtual NM_Status solve(SparseMtrx &A, FloatArray &b, FloatArray &x)=0
Solves the given sparse linear system of equations .
virtual IRResultType initializeFrom(InputRecord *ir)
Initializes receiver according to object description in input reader.
Definition: linearstatic.C:94
LoadBalancer * createLoadBalancer(const char *name, Domain *d)
Definition: classfactory.C:457
std::unique_ptr< TimeStep > currentStep
Current time step.
Definition: engngm.h:231
virtual LoadBalancer * giveLoadBalancer()
Returns reference to receiver&#39;s load balancer.
InputRecord * giveAttributesRecord()
Returns e-model attributes.
Definition: metastep.h:95
virtual TimeStep * giveNextStep()
Returns next time step (next to current step) of receiver.
The reference load vector is obtained as incremental load vector at given time.
Class and object Domain.
Definition: domain.h:115
virtual int giveNumberOfDomainEquations(int di, const UnknownNumberingScheme &num)
Returns number of equations for given domain in active (current time step) time step.
Definition: engngm.C:391
LoadBalancerMonitor * lbm
Definition: engngm.h:293
int giveNumberOfDofManagers() const
Returns number of dof managers in domain.
Definition: domain.h:432
Base class for all matrices stored in sparse format.
Definition: sparsemtrx.h:60
Class representing meta step.
Definition: metastep.h:62
This class implements Newton-Raphson Method, derived from abstract NumericalMethod class for solving ...
Definition: nrsolver.h:95
virtual void solveYourself()
Starts solution process.
SparseMtrx * createSparseMtrx(SparseMtrxType type)
Creates new instance of sparse matrix corresponding to given keyword.
Definition: classfactory.C:105
The purpose of DataStream abstract class is to allow to store/restore context to different streams...
Definition: datastream.h:54
contextIOResultType storeYourself(DataStream &stream) const
Definition: floatarray.C:872
virtual void doStepOutput(TimeStep *tStep)
Prints the ouput of the solution step (using virtual this->printOutputAtservice) to the stream detemi...
Definition: engngm.C:670
std::unique_ptr< TimeStep > previousStep
Previous time step.
Definition: engngm.h:233
double & at(int i)
Coefficient access function.
Definition: floatarray.h:131
The secant stiffness is used and updated only at the beginning of new load step.
Definition: nlinearstatic.h:63
int max(int i, int j)
Returns bigger value form two given decimals.
Definition: mathfem.h:71
MetaStep * giveMetaStep(int i)
Returns the i-th meta step.
Definition: engngm.C:1765
ValueModeType
Type representing the mode of UnknownType or CharType, or similar types.
Definition: valuemodetype.h:78
int nonlocalExt
Flag indicating if nonlocal extension active, which will cause data to be sent between shared element...
Definition: engngm.h:280
The tangent stiffness is used and updated whenever requested.
Definition: nlinearstatic.h:60
long StateCounterType
StateCounterType type used to indicate solution state.
virtual void showSparseMtrxStructure(int type, oofegGraphicContext &gc, TimeStep *tStep)
Shows the sparse structure of required matrix, type == 1 stiffness.
NonLinearStatic_controlType controlMode
Characterizes the type of control used.
#define CM_State
Definition: contextmode.h:46
This base class is an abstraction for all numerical methods solving sparse linear system of equations...
This base class is an abstraction for numerical algorithm.
Definition: nummet.h:80
oofem::oofegGraphicContext gc[OOFEG_LAST_LAYER]
std::unique_ptr< TimeStep > stepWhenIcApply
Solution step when IC (initial conditions) apply.
Definition: engngm.h:229
NonLinearStatic_stiffnessMode stiffMode
EngngModel * giveEngngModel()
Returns engineering model to which receiver is associated.
Definition: domain.C:433
Implementation for assembling tangent matrices in standard monolithic FE-problems.
virtual int estimateMaxPackSize(IntArray &commMap, DataStream &buff, int packUnpackType)
Determines the space necessary for send/receive buffer.
virtual bool hasField(InputFieldType id)=0
Returns true if record contains field identified by idString keyword.
double giveTargetTime()
Returns target time.
Definition: timestep.h:146
bool isParallel() const
Returns true if receiver in parallel mode.
Definition: engngm.h:1056
virtual NumericalMethod * giveNumericalMethod(MetaStep *mStep)
Returns reference to receiver&#39;s numerical method.
Definition: nlinearstatic.C:94
#define _IFT_NonLinearStatic_updateElasticStiffnessFlag
Definition: nlinearstatic.h:54
Base class for dof managers.
Definition: dofmanager.h:113
int giveNumberOfProcesses() const
Returns the number of collaborating processes.
Definition: engngm.h:1060
General IO error.
virtual void solveYourself()
Starts solution process.
Definition: engngm.C:501
#define _IFT_NonLinearStatic_deltat
Definition: nlinearstatic.h:45
StateCounterType internalVarUpdateStamp
Contains last time stamp of internal variable update.
int giveNumber()
Returns domain number.
Definition: domain.h:266
virtual int __giveEquationNumber() const =0
Returns equation number of receiver, usually assigned by emodel.
virtual double giveUnknownComponent(ValueModeType type, TimeStep *tStep, Domain *d, Dof *dof)
Returns requested unknown.
REGISTER_EngngModel(ProblemSequence)
virtual void terminate(TimeStep *tStep)
Terminates the solution of time step.
#define OOFEM_LOG_DEBUG(...)
Definition: logger.h:128
Class implementing an array of integers.
Definition: intarray.h:61
virtual int read(int *data, int count)=0
Reads count integer values into array pointed by data.
virtual void giveInternalForces(FloatArray &answer, bool normFlag, int di, TimeStep *tStep)
Evaluates the nodal representation of internal forces by assembling contributions from individual ele...
FloatArray incrementalLoadVector
Incremental load vector which is applied in loading step (as a whole for direct control or proportion...
#define THROW_CIOERR(e)
Definition: contextioerr.h:61
void printReactionForces(TimeStep *tStep, int id, FILE *out)
Computes and prints reaction forces, computed from nodal internal forces.
NonLinearStatic_controlType
Type determining type of loading control. This type determines the solver to be used.
Definition: nlinearstatic.h:66
NonLinearStatic_stiffnessMode
Type determining the stiffness mode.
Definition: nlinearstatic.h:59
virtual void solveYourselfAt(TimeStep *tStep)
Solves problem for given time step.
#define OOFEM_LOG_RELEVANT(...)
Definition: logger.h:126
void saveStepContext(TimeStep *tStep, ContextMode mode)
Saves context of given solution step, if required (determined using this->giveContextOutputMode() met...
Definition: engngm.C:682
FloatArray initialLoadVectorOfPrescribed
A load vector which does not scale for prescribed DOFs.
FloatArray incrementalLoadVectorOfPrescribed
Incremental Load Vector for prescribed DOFs.
virtual void printStatistics() const
Prints the receiver statistics (one-line) to stdout.
Definition: sparsemtrx.h:275
virtual IRResultType initializeFrom(InputRecord *ir)
Initializes receiver according to object description in input reader.
SparseNonLinearSystemNM::referenceLoadInputModeType refLoadInputMode
The following parameter allows to specify how the reference load vector is obtained from given totalL...
int giveLastStepNumber()
Returns last step number.
Definition: metastep.h:111
double giveTimeIncrement()
Returns solution step associated time increment.
Definition: timestep.h:150
#define _IFT_NonLinearStatic_stiffmode
Definition: nlinearstatic.h:46
void assembleIncrementalReferenceLoadVectors(FloatArray &_incrementalLoadVector, FloatArray &_incrementalLoadVectorOfPrescribed, SparseNonLinearSystemNM::referenceLoadInputModeType _refMode, Domain *sourceDomain, TimeStep *tStep)
virtual int write(const int *data, int count)=0
Writes count integer values from array pointed by data.
virtual void assemble(SparseMtrx &answer, TimeStep *tStep, const MatrixAssembler &ma, const UnknownNumberingScheme &s, Domain *domain)
Assembles characteristic matrix of required type into given sparse matrix.
Definition: engngm.C:803
Abstract base class for all "structural" finite elements.
bool loadBalancingFlag
If set to true, load balancing is active.
Definition: engngm.h:295
MetaStep * giveCurrentMetaStep()
Returns current meta step.
Definition: engngm.C:1684
ProblemCommunicator * nonlocCommunicator
NonLocal Communicator. Necessary when nonlocal constitutive models are used.
Definition: engngm.h:306
int giveNumber()
Returns receiver&#39;s number.
Definition: timestep.h:129
#define OOFEM_LOG_INFO(...)
Definition: logger.h:127
Element * giveElement(int n)
Service for accessing particular domain fe element.
Definition: domain.C:160
void proceedStep(int di, TimeStep *tStep)
NumericalCmpn
Type representing numerical component.
Definition: numericalcmpn.h:46
SparseNonLinearSystemNM * nMethod
Numerical method used to solve the problem.
void doElementOutput(FILE *, TimeStep *)
Does the element output.
#define OOFEM_ERROR(...)
Definition: error.h:61
Callback class for assembling specific types of matrices.
virtual void packMigratingData(TimeStep *tStep)
Packs receiver data when rebalancing load.
int updateSharedDofManagers(FloatArray &answer, const UnknownNumberingScheme &s, int ExchangeTag)
Exchanges necessary remote DofManagers data.
Definition: engngm.C:1957
The reference incremental load vector is defined as totalLoadVector assembled at given time...
Implementation of sparse nonlinear solver with indirect control.
Definition: calmls.h:144
Implementation for assembling external forces vectors in standard monolithic FE-problems.
virtual contextIOResultType saveContext(DataStream &stream, ContextMode mode)
Stores the state of model to output stream.
Abstract base class allowing to control the way, how equations are assigned to individual DOFs...
#define _IFT_EngngModel_initialGuess
Definition: engngm.h:78
contextIOResultType restoreYourself(DataStream &stream)
Definition: floatarray.C:895
#define _IFT_NonLinearStatic_nonlocalext
Definition: nlinearstatic.h:51
ProblemCommunicator * communicator
Communicator.
Definition: engngm.h:303
virtual void printState(FILE *outputStream)
Prints status message of receiver to output stream.
InitialGuess
Means to choose methods for finding a good initial guess.
Definition: engngm.h:197
FloatArray totalDisplacement
Definition: nlinearstatic.h:92
EngngModel * master
Master e-model; if defined receiver is in maintained (slave) mode.
Definition: engngm.h:262
CommunicatorBuff * commBuff
Common Communicator buffer.
Definition: engngm.h:301
int nMetaSteps
Number of meta steps.
Definition: engngm.h:225
virtual int giveNumberOfFirstStep(bool force=false)
Returns number of first time step used by receiver.
Definition: engngm.h:730
virtual void printOutputAt(FILE *file, TimeStep *tStep)
Prints output of receiver to output domain stream, for given time step.
referenceLoadInputModeType
The following parameter allows to specify how the reference load vector is obtained from given totalL...
SparseMtrxType sparseMtrxType
Definition: linearstatic.h:71
bool isEmpty() const
Returns true if receiver is empty.
Definition: floatarray.h:222
The initial elastic stiffness is used in the whole solution.
Definition: nlinearstatic.h:62
double giveIntrinsicTime()
Returns intrinsic time, e.g. time in which constitutive model is evaluated.
Definition: timestep.h:148
virtual void updateDomainLinks()
Updates domain links after the domains of receiver have changed.
Definition: linearstatic.C:300
#define _IFT_NonLinearStatic_keepll
Definition: nlinearstatic.h:48
double deltaT
Intrinsic time increment.
Abstract base class representing general load balancer.
Definition: loadbalancer.h:108
Describes the direct control where load or displacement (or both) are controlled. ...
Definition: nlinearstatic.h:68
virtual contextIOResultType saveContext(DataStream &stream, ContextMode mode)
Stores the state of model to output stream.
Definition: engngm.C:1527
#define NM_None
Definition: nmstatus.h:46
Initializes the variable VERBOSE, in order to get a few intermediate messages on screen: beginning an...
virtual void updateComponent(TimeStep *tStep, NumericalCmpn, Domain *d)
Updates components mapped to numerical method if necessary during solution process.
Class representing vector of real numbers.
Definition: floatarray.h:82
int estimatePackSize(DataStream &buff)
Estimates the necessary pack size to hold all packed data of receiver.
Definition: element.C:1575
#define _IFT_NonLinearStatic_refloadmode
Definition: nlinearstatic.h:47
virtual void updateDomainLinks()
Updates domain links after the domains of receiver have changed.
virtual void reinitialize()
virtual void computeExternalLoadReactionContribution(FloatArray &reactions, TimeStep *tStep, int di)
Computes the contribution external loading to reaction forces in given domain.
FloatArray incrementOfDisplacement
Definition: nlinearstatic.h:92
IRResultType
Type defining the return values of InputRecord reading operations.
Definition: irresulttype.h:47
virtual int givePackSizeOfDouble(int count)=0
std::unique_ptr< SparseMtrx > stiffnessMatrix
Definition: linearstatic.h:66
No special treatment for new iterations. Probably means ending up using for all free dofs...
Definition: engngm.h:198
int giveMetaStepNumber()
Returns receiver&#39;s meta step number.
Definition: timestep.h:135
#define _IFT_NonLinearStatic_controlmode
Definition: nlinearstatic.h:44
int giveNumberOfTimeStepWhenIcApply()
Returns the time step number, when initial conditions should apply.
Definition: engngm.h:754
InitialGuess initialGuessType
The initial guess type to use before starting the nonlinear solver.
CharType
Definition: chartype.h:87
NonLinearStatic(int i, EngngModel *_master=NULL)
Definition: nlinearstatic.C:66
virtual void reinitialize()
Reinitializes the receiver.
Definition: nummet.h:112
Class representing the general Input Record.
Definition: inputrecord.h:101
The representation of EngngModel default prescribed unknown numbering.
#define _IFT_NonLinearStatic_nonlocstiff
Definition: nlinearstatic.h:50
virtual contextIOResultType restoreContext(DataStream &stream, ContextMode mode)
Restores the state of model from output stream.
void zero()
Zeroes all coefficients of receiver.
Definition: floatarray.C:658
This class implements linear static engineering problem.
Definition: linearstatic.h:63
Class implementing single timer, providing wall clock and user time capabilities. ...
Definition: timer.h:46
void assemblePrescribedExtrapolatedForces(FloatArray &answer, TimeStep *tStep, CharType type, Domain *domain)
Definition: engngm.C:1417
int giveRank() const
Returns domain rank in a group of collaborating processes (0..groupSize-1)
Definition: engngm.h:1058
int giveNumber()
Returns receiver&#39;s number.
Definition: metastep.h:89
void initializeCommMaps(bool forceInit=false)
Definition: engngm.C:1942
ClassFactory & classFactory
Definition: classfactory.C:59
long ContextMode
Context mode (mask), defining the type of information written/read to/from context.
Definition: contextmode.h:43
LoadBalancerMonitor * createLoadBalancerMonitor(const char *name, EngngModel *e)
Definition: classfactory.C:447
LoadBalancer * lb
Load Balancer.
Definition: engngm.h:292
#define CM_Definition
Definition: contextmode.h:47
virtual ErrorEstimator * giveDomainErrorEstimator(int n)
Service for accessing ErrorEstimator corresponding to particular domain.
Definition: engngm.h:349
OutputManager * giveOutputManager()
Returns domain output manager.
Definition: domain.C:1436
int giveVersion()
Returns receiver&#39;s version.
Definition: timestep.h:133
virtual TimeStep * giveSolutionStepWhenIcApply(bool force=false)
Returns the solution step when Initial Conditions (IC) apply.
Definition: engngm.h:720
virtual void assemble(SparseMtrx &answer, TimeStep *tStep, const MatrixAssembler &ma, const UnknownNumberingScheme &, Domain *domain)
Assembles characteristic matrix of required type into given sparse matrix.
virtual LoadBalancerMonitor * giveLoadBalancerMonitor()
Returns reference to receiver&#39;s load balancer monitor.
FloatArray internalForcesEBENorm
Norm of nodal internal forces evaluated on element by element basis (squared)
A generalized norm of displacement and loading vectors is controlled. In current implementation, the CALM solver is used, the reference load vector is FIXED.
Definition: nlinearstatic.h:67
The Communicator and corresponding buffers (represented by this class) are separated in order to allo...
Definition: communicator.h:60
std::vector< std::unique_ptr< Element > > & giveElements()
Definition: domain.h:279
virtual void setDomain(Domain *d)
sets associated Domain
Definition: loadbalancer.h:162
Abstract base class representing the "problem" under consideration.
Definition: engngm.h:181
virtual void setDomain(Domain *d)
Definition: nummet.h:114
#define IR_GIVE_OPTIONAL_FIELD(__ir, __value, __id)
Macro facilitating the use of input record reading methods.
Definition: inputrecord.h:78
virtual TimeStep * giveCurrentStep(bool force=false)
Returns current time step.
Definition: engngm.h:683
The secant stiffness is used and updated whenever requested.
Definition: nlinearstatic.h:61
virtual contextIOResultType restoreContext(DataStream &stream, ContextMode mode)
Restores the state of model from output stream.
Definition: engngm.C:1592
virtual void updateAttributes(MetaStep *mStep)
Update receiver attributes according to step metaStep attributes.
the oofem namespace is to define a context or scope in which all oofem names are defined.
bool isNotEmpty() const
Returns true if receiver is not empty.
Definition: floatarray.h:220
virtual void updateLoadVectors(TimeStep *tStep)
Domain * giveDomain(int n)
Service for accessing particular problem domain.
Definition: engngm.C:1720
problemScale pScale
Multiscale mode.
Definition: engngm.h:257
FloatArray initialLoadVector
A load vector already applied, which does not scales.
Definition: nlinearstatic.h:96
Abstract class Dof represents Degree Of Freedom in finite element mesh.
Definition: dof.h:93
DofManager * giveDofManager(int n)
Service for accessing particular domain dof manager.
Definition: domain.C:314
void negated()
Switches the sign of every coefficient of receiver.
Definition: floatarray.C:739
void doDofManOutput(FILE *, TimeStep *)
Does the dofmanager output.
Definition: outputmanager.C:78
void startTimer()
Definition: timer.C:69
double getUtime()
Returns total user time elapsed in seconds.
Definition: timer.C:105
virtual TimeStep * giveSolutionStepWhenIcApply(bool force=false)
Returns the solution step when Initial Conditions (IC) apply.
void assembleVector(FloatArray &answer, TimeStep *tStep, const VectorAssembler &va, ValueModeType mode, const UnknownNumberingScheme &s, Domain *domain, FloatArray *eNorms=NULL)
Assembles characteristic vector of required type from dofManagers, element, and active boundary condi...
Definition: engngm.C:986
Class representing solution step.
Definition: timestep.h:80
DofManager is shared by neighboring partitions, it is necessary to sum contributions from all contrib...
Definition: dofmanager.h:82
void stopTimer()
Definition: timer.C:77
void add(const FloatArray &src)
Adds array src to receiver.
Definition: floatarray.C:156
Abstract base class representing general load balancer monitor.
Definition: loadbalancer.h:68
virtual NM_Status solve(SparseMtrx &K, FloatArray &R, FloatArray *R0, FloatArray &X, FloatArray &dX, FloatArray &F, const FloatArray &internalForcesEBENorm, double &s, referenceLoadInputModeType rlm, int &nite, TimeStep *tStep)=0
Solves the given sparse linear system of equations .
virtual void unpackMigratingData(TimeStep *tStep)
Unpacks receiver data when rebalancing load.
#define _IFT_NonLinearStatic_donotfixload
Definition: nlinearstatic.h:49
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:30 for OOFEM by doxygen 1.8.11 written by Dimitri van Heesch, © 1997-2011