MAST
Multidisciplinary-design Adaptation and Sensitivity Toolkit (MAST)
transient_solver_base.cpp
Go to the documentation of this file.
1 /*
2  * MAST: Multidisciplinary-design Adaptation and Sensitivity Toolkit
3  * Copyright (C) 2013-2020 Manav Bhatia and MAST authors
4  *
5  * This library is free software; you can redistribute it and/or
6  * modify it under the terms of the GNU Lesser General Public
7  * License as published by the Free Software Foundation; either
8  * version 2.1 of the License, or (at your option) any later version.
9  *
10  * This library is distributed in the hope that it will be useful,
11  * but WITHOUT ANY WARRANTY; without even the implied warranty of
12  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
13  * Lesser General Public License for more details.
14  *
15  * You should have received a copy of the GNU Lesser General Public
16  * License along with this library; if not, write to the Free Software
17  * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
18  */
19 
20 // MAST includes
23 #include "base/assembly_base.h"
25 #include "base/nonlinear_system.h"
27 
28 // libMesh includes
29 #include "libmesh/numeric_vector.h"
30 #include "libmesh/dof_map.h"
31 #include "libmesh/sparse_matrix.h"
32 #include "libmesh/linear_solver.h"
33 
34 
35 
37  unsigned int n):
39 dt (0.),
40 _first_step (true),
41 _first_sensitivity_step (true),
42 _ode_order (o),
43 _n_iters_to_store (n),
44 _assembly_ops (nullptr),
45 _if_highest_derivative_solution (false) {
46 
47 }
48 
49 
50 
51 
53 
55 }
56 
57 
58 
59 void
61 
62  libmesh_assert(_assembly_ops);
64  _assembly_ops->set_assembly(assembly);
65 }
66 
67 
68 void
70 
72  if (_assembly_ops)
73  _assembly_ops->clear_assembly();
74 }
75 
76 
79 
80  libmesh_assert(_assembly_ops);
81 
82  return *_assembly_ops;
83 }
84 
85 
86 void
89 
90  libmesh_assert(_system);
91  libmesh_assert(_discipline);
92 
94  &sys = _system->system();
95 
96  // make sure that the previous assembly association has been cleared
98 
99  _assembly_ops = &assembly_ops;
100 
101  // now, add the vectors
102  std::string nm;
103  for (unsigned int i=0; i<_n_iters_to_store; i++) {
104  std::ostringstream iter;
105  iter << i;
106 
107  // add the solution only for previous iterations. The current
108  // solution is obtained from system->solution
109  if (i>0) {
110 
111  nm = "transient_solution_";
112  nm += iter.str();
113  sys.add_vector(nm, true, libMesh::GHOSTED);
114  }
115 
116  nm = "transient_solution_sensitivity_";
117  nm += iter.str();
118  sys.add_vector(nm, true, libMesh::GHOSTED);
119 
120  // add the velocity
121  nm = "transient_velocity_";
122  nm += iter.str();
123  sys.add_vector(nm, true, libMesh::GHOSTED);
124  nm = "transient_velocity_sensitivity_";
125  nm += iter.str();
126  sys.add_vector(nm, true, libMesh::GHOSTED);
127 
128  if (_ode_order > 1) {
129  // add the acceleration
130  nm = "transient_acceleration_";
131  nm += iter.str();
132  sys.add_vector(nm, true, libMesh::GHOSTED);
133  nm = "transient_acceleration_sensitivity_";
134  nm += iter.str();
135  sys.add_vector(nm, true, libMesh::GHOSTED);
136  }
137  }
138 
139  _first_step = true;
140 }
141 
142 
143 
144 void
146 
147  // if no system has been set so far, nothing needs to be done
148  if (!_system)
149  return;
150 
152  &sys = _system->system();
153 
154  // clear the transient solutions stored in system for solution
155  // number of time steps to store
156  std::string nm;
157  for (unsigned int i=0; i<_n_iters_to_store; i++) {
158  std::ostringstream iter;
159  iter << i;
160 
161  // remove the solution
162  nm = "transient_solution_";
163  nm += iter.str();
164  if (sys.have_vector(nm)) sys.remove_vector(nm);
165  nm = "transient_solution_sensitivity_";
166  nm += iter.str();
167  if (sys.have_vector(nm)) sys.remove_vector(nm);
168 
169  // remove the velocity
170  nm = "transient_velocity_";
171  nm += iter.str();
172  if (sys.have_vector(nm)) sys.remove_vector(nm);
173  nm = "transient_velocity_sensitivity_";
174  nm += iter.str();
175  if (sys.have_vector(nm)) sys.remove_vector(nm);
176 
177  if (_ode_order > 1) {
178  // remove the acceleration
179  nm = "transient_acceleration_";
180  nm += iter.str();
181  if (sys.have_vector(nm)) sys.remove_vector(nm);
182  nm = "transient_acceleration_sensitivity_";
183  nm += iter.str();
184  if (sys.have_vector(nm)) sys.remove_vector(nm);
185 
186  }
187  }
188 
189  _assembly_ops = nullptr;
190  _first_step = true;
191 }
192 
193 
194 
195 
196 libMesh::NumericVector<Real>&
197 MAST::TransientSolverBase::solution(unsigned int prev_iter) const {
198 
199  // make sura that prev_iter is within acceptable bounds
200  libmesh_assert_less(prev_iter, _n_iters_to_store);
201 
202  // make sure that the vectors have been initialized
203  std::ostringstream oss;
204  oss << prev_iter;
205 
206  if (prev_iter) {
207  // get references to the solution
208  std::string
209  nm = "transient_solution_" + oss.str();
210 
211  return _system->system().get_vector(nm);
212  }
213  else
214  return *_system->system().current_local_solution;
215 }
216 
217 
218 
219 libMesh::NumericVector<Real>&
220 MAST::TransientSolverBase::solution_sensitivity(unsigned int prev_iter) const {
221 
222  // make sura that prev_iter is within acceptable bounds
223  libmesh_assert_less(prev_iter, _n_iters_to_store);
224 
225  // make sure that the vectors have been initialized
226  std::ostringstream oss;
227  oss << prev_iter;
228 
229  // get references to the solution
230  std::string
231  nm = "transient_solution_sensitivity_" + oss.str();
232 
233  return _system->system().get_vector(nm);
234 }
235 
236 
237 
238 libMesh::NumericVector<Real>&
239 MAST::TransientSolverBase::velocity(unsigned int prev_iter) const {
240 
241  // make sura that prev_iter is within acceptable bounds
242  libmesh_assert_less(prev_iter, _n_iters_to_store);
243 
244  // make sure that the vectors have been initialized
245  std::ostringstream oss;
246  oss << prev_iter;
247 
248  // get references to the solution
249  std::string
250  nm = "transient_velocity_" + oss.str();
251 
252  return _system->system().get_vector(nm);
253 }
254 
255 
256 
257 libMesh::NumericVector<Real>&
258 MAST::TransientSolverBase::velocity_sensitivity(unsigned int prev_iter) const {
259 
260  // make sura that prev_iter is within acceptable bounds
261  libmesh_assert_less(prev_iter, _n_iters_to_store);
262 
263  // make sure that the vectors have been initialized
264  std::ostringstream oss;
265  oss << prev_iter;
266 
267  // get references to the solution
268  std::string
269  nm = "transient_velocity_sensitivity_" + oss.str();
270 
271  return _system->system().get_vector(nm);
272 }
273 
274 
275 
276 
277 libMesh::NumericVector<Real>&
278 MAST::TransientSolverBase::acceleration(unsigned int prev_iter) const {
279 
280  // make sura that prev_iter is within acceptable bounds
281  libmesh_assert_less(prev_iter, _n_iters_to_store);
282 
283  // make sure that the vectors have been initialized
284  std::ostringstream oss;
285  oss << prev_iter;
286 
287  // get references to the solution
288  std::string
289  nm = "transient_acceleration_" + oss.str();
290 
291  return _system->system().get_vector(nm);
292 }
293 
294 
295 
296 libMesh::NumericVector<Real>&
298 
299  // make sura that prev_iter is within acceptable bounds
300  libmesh_assert_less(prev_iter, _n_iters_to_store);
301 
302  // make sure that the vectors have been initialized
303  std::ostringstream oss;
304  oss << prev_iter;
305 
306  // get references to the solution
307  std::string
308  nm = "transient_acceleration_sensitivity_" + oss.str();
309 
310  return _system->system().get_vector(nm);
311 }
312 
313 
314 
315 void
317 
318  // make sure that the system has been specified
319  libmesh_assert_msg(_system, "System pointer is nullptr.");
320 
321  // ask the Newton solver to solve for the system solution
322  _system->system().solve(*this, assembly);
323 
324 }
325 
326 
327 
328 void
330  const MAST::FunctionBase& f) {
331 
332  // make sure that the system has been specified
333  libmesh_assert_msg(_system, "System pointer is nullptr.");
334 
335  // ask the Newton solver to solve for the system solution
336  _system->system().sensitivity_solve(this->solution(), false, *this, assembly, f);
337 
338  this->solution_sensitivity() = assembly.system().get_sensitivity_solution();
339 }
340 
341 
342 
343 void
346  bool increment_time) {
347 
348  libmesh_assert(_first_step);
349  libmesh_assert(_system);
350  libmesh_assert(_discipline);
351  libmesh_assert(!_assembly);
352 
353  // tell the solver that the current solution being obtained is for the
354  // highest time derivative
355  _if_highest_derivative_solution = true;
356 
358  &sys = _system->system();
359 
360  // Build the residual and Jacobian
361  assembly.set_elem_operation_object(*this);
362  assembly.residual_and_jacobian(*sys.solution, sys.rhs, sys.matrix, sys);
363  assembly.clear_elem_operation_object();
364 
365  // reset the solution flag
366  _if_highest_derivative_solution = false;
367 
368  std::pair<unsigned int, Real>
369  solver_params = sys.get_linear_solve_parameters();
370 
371  // Solve the linear system.
372  libMesh::SparseMatrix<Real> *
373  pc = sys.request_matrix("Preconditioner");
374 
375  std::unique_ptr<libMesh::NumericVector<Real> >
376  dvec(velocity().zero_clone().release());
377 
378  sys.linear_solver->solve (*sys.matrix, pc,
379  *dvec,
380  *sys.rhs,
381  solver_params.second,
382  solver_params.first);
383 
384  libMesh::NumericVector<Real> *vec = nullptr;
385 
386  switch (_ode_order) {
387 
388  case 1:
389  vec = &velocity();
390  break;
391 
392  case 2:
393  vec = &acceleration();
394  break;
395 
396  default:
397  // higher than 2 derivative not implemented yet.
398  libmesh_error();
399  }
400 
401  vec->add(-1., *dvec);
402 
403  // The linear solver may not have fit our constraints exactly
404 #ifdef LIBMESH_ENABLE_CONSTRAINTS
405  sys.get_dof_map().enforce_constraints_exactly(sys, vec, /* homogeneous = */ true);
406 #endif
407 
408  sys.update();
409 
410  // next, move all the solutions and velocities into older
411  // time step locations
412  for (unsigned int i=_n_iters_to_store-1; i>0; i--) {
413  this->solution(i) = this->solution(i-1);
414  this->velocity(i) = this->velocity(i-1);
415 
416  if (_ode_order > 1)
417  this->acceleration(i) = this->acceleration(i-1);
418  }
419 
420  // finally, update the system time
421  if (increment_time) sys.time += dt;
422  _first_step = false;
423 }
424 
425 
426 
427 
428 void
431  (MAST::AssemblyBase& assembly,
432  const MAST::FunctionBase& f) {
433 
434  libmesh_assert(_first_sensitivity_step);
435  libmesh_assert(_system);
436  libmesh_assert(_discipline);
437  libmesh_assert(!_assembly);
438 
439  // tell the solver that the current solution being obtained is for the
440  // highest time derivative
441  _if_highest_derivative_solution = true;
442 
444  &sys = _system->system();
445 
446  // Build the Jacobian
447  assembly.set_elem_operation_object(*this);
448  assembly.residual_and_jacobian(*sys.solution, nullptr, sys.matrix, sys);
449  // this should assembl the sensitivity of the residual
450  assembly.sensitivity_assemble(*sys.solution, true, f, *sys.rhs);
451  assembly.clear_elem_operation_object();
452 
453  // reset the solution flag
454  _if_highest_derivative_solution = false;
455 
456  std::pair<unsigned int, Real>
457  solver_params = sys.get_linear_solve_parameters();
458 
459  // Solve the linear system.
460  libMesh::SparseMatrix<Real> *
461  pc = sys.request_matrix("Preconditioner");
462 
463  std::unique_ptr<libMesh::NumericVector<Real> >
464  dvec(velocity().zero_clone().release());
465 
466  sys.linear_solver->solve (*sys.matrix, pc,
467  *dvec,
468  *sys.rhs,
469  solver_params.second,
470  solver_params.first);
471 
472  libMesh::NumericVector<Real> *vec = nullptr;
473 
474  switch (_ode_order) {
475 
476  case 1:
477  vec = &velocity_sensitivity();
478  break;
479 
480  case 2:
481  vec = &acceleration_sensitivity();
482  break;
483 
484  default:
485  // higher than 2 derivative not implemented yet.
486  libmesh_error();
487  }
488 
489  vec->add(-1., *dvec);
490 
491  // The linear solver may not have fit our constraints exactly
492 #ifdef LIBMESH_ENABLE_CONSTRAINTS
493  sys.get_dof_map().enforce_constraints_exactly(sys, vec, /* homogeneous = */ true);
494 #endif
495 
496  // next, move all the solutions and velocities into older
497  // time step locations
498  for (unsigned int i=_n_iters_to_store-1; i>0; i--) {
499 
500  this->solution_sensitivity(i) = this->solution_sensitivity(i-1);
501  this->velocity_sensitivity(i) = this->velocity_sensitivity(i-1);
502 
503  if (_ode_order > 1)
505  }
506 
507  // finally, update the system time
508  sys.time += dt;
509  _first_sensitivity_step = false;
510 }
511 
512 
513 
514 void
516 build_local_quantities(const libMesh::NumericVector<Real>& current_sol,
517  std::vector<libMesh::NumericVector<Real>*>& sol) {
518 
519  // make sure that the system has been specified
520  libmesh_assert(_system);
521 
523  &sys = _system->system();
524 
525  // make sure there are no solutions in sol
526  libmesh_assert(!sol.size());
527 
528  // resize the vector to store the quantities
529  sol.resize(_ode_order+1);
530 
531  const std::vector<libMesh::dof_id_type>&
532  send_list = sys.get_dof_map().get_send_list();
533 
534 
535  for ( unsigned int i=0; i<=_ode_order; i++) {
536 
537  sol[i] = libMesh::NumericVector<Real>::build(sys.comm()).release();
538  sol[i]->init(sys.n_dofs(),
539  sys.n_local_dofs(),
540  send_list,
541  false,
542  libMesh::GHOSTED);
543 
544  switch (i) {
545 
546  case 0: {
547 
548  if (!_if_highest_derivative_solution)
549  // copy the solution to this
550  current_sol.localize(*sol[i], send_list);
551  else
552  solution().localize(*sol[i], send_list);
553  }
554  break;
555 
556  case 1: {
557 
558  // update the current local velocity vector
559  libMesh::NumericVector<Real>& vel = this->velocity();
560 
561  if (!_if_highest_derivative_solution)
562  // calculate the velocity and localize it
563  // we use the localized sol since libMesh vector
564  // algebra needs it for ghosted vectors.
565  update_velocity(vel, *sol[0]);
566 
567  vel.localize(*sol[i], send_list);
568  }
569  break;
570 
571  case 2: {
572 
573  // update the current local acceleration vector
574  libMesh::NumericVector<Real>& acc = this->acceleration();
575 
576  if (!_if_highest_derivative_solution)
577  // calculate the acceleration and localize it
578  // we use the localized sol since libMesh vector
579  // algebra needs it for ghosted vectors.
580  update_acceleration(acc, *sol[0]);
581 
582  acc.localize(*sol[i], send_list);
583  }
584  break;
585 
586  default:
587  // should not get here
588  libmesh_error();
589  break;
590  }
591  }
592 }
593 
594 
595 
596 void
598 build_sensitivity_local_quantities(unsigned int prev_iter,
599  std::vector<libMesh::NumericVector<Real>*>& sol) {
600 
601  // make sure that the system has been specified
602  libmesh_assert(_system);
603  libmesh_assert_less_equal(prev_iter, _n_iters_to_store);
604 
606  &sys = _system->system();
607 
608  // make sure there are no solutions in sol
609  libmesh_assert(!sol.size());
610 
611  // resize the vector to store the quantities
612  sol.resize(_ode_order+1);
613 
614  const std::vector<libMesh::dof_id_type>&
615  send_list = sys.get_dof_map().get_send_list();
616 
617 
618  for ( unsigned int i=0; i<=_ode_order; i++) {
619 
620  sol[i] = libMesh::NumericVector<Real>::build(sys.comm()).release();
621  sol[i]->init(sys.n_dofs(),
622  sys.n_local_dofs(),
623  send_list,
624  false,
625  libMesh::GHOSTED);
626 
627  switch (i) {
628 
629  case 0:
630  solution_sensitivity(prev_iter).localize(*sol[i], send_list);
631  break;
632 
633  case 1:
634  velocity_sensitivity(prev_iter).localize(*sol[i], send_list);
635  break;
636 
637  case 2:
638  acceleration_sensitivity(prev_iter).localize(*sol[i], send_list);
639  break;
640 
641  default:
642  // should not get here
643  libmesh_error();
644  break;
645  }
646  }
647 }
648 
649 
650 
651 
652 
653 void
654 MAST::TransientSolverBase::
655 build_perturbed_local_quantities(const libMesh::NumericVector<Real>& current_dsol,
656  std::vector<libMesh::NumericVector<Real>*>& sol) {
657 
658  // make sure that the system has been specified
659  libmesh_assert(_system);
660 
662  &sys = _system->system();
663 
664  // make sure there are no solutions in sol
665  libmesh_assert(!sol.size());
666 
667  // resize the vector to store the quantities
668  sol.resize(_ode_order+1);
669 
670  const std::vector<libMesh::dof_id_type>&
671  send_list = sys.get_dof_map().get_send_list();
672 
673 
674  std::unique_ptr<libMesh::NumericVector<Real> >
675  tmp(this->solution().zero_clone().release());
676 
677  for ( unsigned int i=0; i<=_ode_order; i++) {
678 
679  sol[i] = libMesh::NumericVector<Real>::build(sys.comm()).release();
680  sol[i]->init(sys.n_dofs(),
681  sys.n_local_dofs(),
682  send_list,
683  false,
684  libMesh::GHOSTED);
685 
686  switch (i) {
687 
688  case 0: {
689 
690  current_dsol.localize(*sol[i], send_list);
691  if (_if_highest_derivative_solution)
692  // only the highest derivative is perturbed since
693  // all lower derivative quantities are known
694  sol[i]->zero();
695  }
696  break;
697 
698  case 1: {
699 
700  if (_ode_order == 1 && _if_highest_derivative_solution)
701  // the provided solution is the current velocity increment
702  current_dsol.localize(*sol[i], send_list);
703  else if (_if_highest_derivative_solution)
704  sol[i]->zero();
705  else {
706  update_delta_velocity(*tmp, current_dsol);
707  tmp->localize(*sol[i], send_list);
708  }
709  }
710  break;
711 
712  case 2: {
713 
714  if (_ode_order == 2 && _if_highest_derivative_solution)
715  // the provided solution is the current velocity increment
716  current_dsol.localize(*sol[i], send_list);
717  else if (_if_highest_derivative_solution)
718  sol[i]->zero();
719  else {
720  update_delta_acceleration(*tmp, current_dsol);
721  tmp->localize(*sol[i], send_list);
722  }
723  }
724  break;
725 
726  default:
727  // should not get here
728  libmesh_error();
729  break;
730  }
731  }
732 }
733 
734 
735 
736 
737 void
739 
740  libmesh_assert(_system);
741  libmesh_assert(_discipline);
742 
744  &sys = _system->system();
745 
746  // first ask the solver to update the velocity and acceleration vector
747  sys.update();
748  update_velocity(this->velocity(), *sys.current_local_solution);
749 
750  if (_ode_order > 1)
751  update_acceleration(this->acceleration(), *sys.current_local_solution);
752 
753  // next, move all the solutions and velocities into older
754  // time step locations
755  for (unsigned int i=_n_iters_to_store-1; i>0; i--) {
756  this->solution(i) = this->solution(i-1);
757  this->velocity(i) = this->velocity(i-1);
758 
759  if (_ode_order > 1)
760  this->acceleration(i) = this->acceleration(i-1);
761  }
762 
763  // finally, update the system time
764  if (increment_time) sys.time += dt;
765  _first_step = false;
766 }
767 
768 
769 
770 void
772 
774  &sys = _system->system();
775 
776  // first ask the solver to update the velocity and acceleration vector
778 
779  if (_ode_order > 1)
781 
782  // next, move all the solutions and velocities into older
783  // time step locations
784  for (unsigned int i=_n_iters_to_store-1; i>0; i--) {
785  this->solution_sensitivity(i) = this->solution_sensitivity(i-1);
786  this->velocity_sensitivity(i) = this->velocity_sensitivity(i-1);
787 
788  if (_ode_order > 1)
790  }
791 
792  // finally, update the system time
793  sys.time += dt;
794  _first_sensitivity_step = false;
795 }
796 
797 
798 
799 void
801  const libMesh::Elem& ref_elem,
802  MAST::GeomElem &elem) const {
803 
804  libmesh_assert(_assembly_ops);
805  _assembly_ops->set_elem_data(dim, ref_elem, elem);
806 }
807 
808 
809 void
811 
812  libmesh_assert(_assembly_ops);
813  _assembly_ops->init(elem);
814 }
815 
816 
817 
818 
819 void
821 
822  libmesh_assert(_assembly_ops);
823  _assembly_ops->clear_elem();
824 }
825 
826 
const MAST::NonlinearSystem & system() const
virtual void build_local_quantities(const libMesh::NumericVector< Real > &current_sol, std::vector< libMesh::NumericVector< Real > *> &qtys)
localizes the relevant solutions for system assembly.
MAST::NonlinearSystem & system()
virtual void advance_time_step_with_sensitivity()
advances the time step and copies the current sensitivity solution to old sensitivity solution...
virtual void advance_time_step(bool increment_time=true)
advances the time step and copies the current solution to old solution, and so on.
virtual void sensitivity_solve(MAST::AssemblyBase &assembly, const MAST::FunctionBase &f)
solvers the current time step for sensitivity wrt f
virtual void update_sensitivity_acceleration(libMesh::NumericVector< Real > &acc, const libMesh::NumericVector< Real > &sol)=0
update the transient sensitivity acceleration based on the current sensitivity solution ...
This class implements a system for solution of nonlinear systems.
virtual void set_assembly(MAST::AssemblyBase &assembly)
sets the assembly object
virtual MAST::TransientAssemblyElemOperations & get_elem_operation_object()
virtual void clear_elem_operation_object()
clears the association of this object with the assembly element operation object. ...
virtual void update_delta_acceleration(libMesh::NumericVector< Real > &acc, const libMesh::NumericVector< Real > &sol)=0
update the perturbation in transient acceleration based on the current perturbed solution ...
virtual void set_assembly(MAST::AssemblyBase &assembly)
sets the assembly object
virtual void update_delta_velocity(libMesh::NumericVector< Real > &vel, const libMesh::NumericVector< Real > &sol)=0
update the perturbation in transient velocity based on the current perturbed solution ...
libMesh::NumericVector< Real > & solution(unsigned int prev_iter=0) const
virtual void set_elem_operation_object(MAST::AssemblyElemOperations &elem_ops)
attaches a element operation to this object, and associated this with the element operation object...
MAST::PhysicsDisciplineBase * _discipline
virtual void set_elem_operation_object(MAST::TransientAssemblyElemOperations &elem_ops)
Attaches the assembly elem operations object that provides the x_dot, M and J quantities for the elem...
libMesh::NumericVector< Real > & solution_sensitivity(unsigned int prev_iter=0) const
libMesh::NumericVector< Real > & acceleration_sensitivity(unsigned int prev_iter=0) const
virtual void set_elem_data(unsigned int dim, const libMesh::Elem &ref_elem, MAST::GeomElem &elem) const =0
some analyses may want to set additional element data before initialization of the GeomElem...
virtual void update_velocity(libMesh::NumericVector< Real > &vel, const libMesh::NumericVector< Real > &sol)=0
update the transient velocity based on the current solution
libMesh::NumericVector< Real > & velocity(unsigned int prev_iter=0) const
virtual void clear_assembly()
clears the assembly object
libMesh::NumericVector< Real > & acceleration(unsigned int prev_iter=0) const
virtual std::pair< unsigned int, Real > get_linear_solve_parameters()
calls NonlinearImplicitSystem::set_solver_parameters() before accessing the values.
virtual void init(const MAST::GeomElem &elem)=0
initializes the object for calculation of element quantities for the specified elem.
virtual void update_acceleration(libMesh::NumericVector< Real > &acc, const libMesh::NumericVector< Real > &sol)=0
update the transient acceleration based on the current solution
void solve_highest_derivative_and_advance_time_step(MAST::AssemblyBase &assembly, bool increment_time=true)
To be used only for initial conditions.
virtual void clear_assembly()
clears the assembly object
This class acts as a wrapper around libMesh::Elem for the purpose of providing a uniform interface fo...
Definition: geom_elem.h:59
virtual void solve(MAST::AssemblyBase &assembly)
solves the current time step for solution and velocity
std::unique_ptr< libMesh::LinearSolver< Real > > linear_solver
The LinearSolver for solution of the linear equations.
virtual void clear_elem()
clears the element initialization
virtual void residual_and_jacobian(const libMesh::NumericVector< Real > &X, libMesh::NumericVector< Real > *R, libMesh::SparseMatrix< Real > *J, libMesh::NonlinearImplicitSystem &S)
function that assembles the matrices and vectors quantities for nonlinear solution ...
virtual void build_sensitivity_local_quantities(unsigned int prev_iter, std::vector< libMesh::NumericVector< Real > *> &qtys)
localizes the relevant solutions for system assembly.
libMesh::NumericVector< Real > & velocity_sensitivity(unsigned int prev_iter=0) const
virtual void sensitivity_solve(const libMesh::NumericVector< Real > &X, bool if_localize_sol, MAST::AssemblyElemOperations &elem_ops, MAST::AssemblyBase &assembly, const MAST::FunctionBase &p, bool if_assemble_jacobian=true)
Solves the sensitivity problem for the provided parameter.
TransientSolverBase(unsigned int o, unsigned int n)
constructor requires the number of iterations to store for the derived solver.
virtual void update_sensitivity_velocity(libMesh::NumericVector< Real > &vel, const libMesh::NumericVector< Real > &sol)=0
update the transient sensitivity velocity based on the current sensitivity solution ...
void solve_highest_derivative_and_advance_time_step_with_sensitivity(MAST::AssemblyBase &assembly, const MAST::FunctionBase &f)
solves for the sensitivity of highest derivative and advances the time-step.
virtual bool sensitivity_assemble(const libMesh::NumericVector< Real > &X, bool if_localize_sol, const MAST::FunctionBase &f, libMesh::NumericVector< Real > &sensitivity_rhs, bool close_vector=true)
Assembly function.
virtual void clear_elem_operation_object()
Clears the assembly elem operations object.
virtual void solve(MAST::AssemblyElemOperations &elem_ops, MAST::AssemblyBase &assembly)
solves the nonlinear problem with the specified assembly operation object
MAST::SystemInitialization * _system