explicit_system.C
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2014 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
3 
4 // This library is free software; you can redistribute it and/or
5 // modify it under the terms of the GNU Lesser General Public
6 // License as published by the Free Software Foundation; either
7 // version 2.1 of the License, or (at your option) any later version.
8 
9 // This library is distributed in the hope that it will be useful,
10 // but WITHOUT ANY WARRANTY; without even the implied warranty of
11 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
12 // Lesser General Public License for more details.
13 
14 // You should have received a copy of the GNU Lesser General Public
15 // License along with this library; if not, write to the Free Software
16 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 
18 
19 
20 // C++ includes
21 
22 // Local includes
24 #include "libmesh/numeric_vector.h"
25 
26 namespace libMesh
27 {
28 
29 
30 // ------------------------------------------------------------
31 // ExplicitSystem implementation
33  const std::string& name_in,
34  const unsigned int number_in) :
35  Parent (es, name_in, number_in),
36  rhs(NULL)
37 
38 {
39  //rhs = &(this->add_vector ("RHS Vector"));
40 }
41 
42 
43 
45 {
46  // clear data
47  this->clear();
48 }
49 
50 
51 
53 {
54  // Clear the parent data
55  Parent::clear();
56 
57  // NULL-out the vector. Note that
58  // System::clear() actually deleted it.
59  rhs = NULL;
60 }
61 
62 
63 
65 {
66  // Add the system RHS.
67  // (We must do this before initializing the System data,
68  // then we lose the opportunity to add vectors).
69  this->add_system_rhs ();
70 
71  // initialize parent data
73 }
74 
75 
76 
78 {
79  // initialize parent data
81 
82  // not necessary, handled by the parent!
83  // Resize the RHS conformal to the current mesh
84  //rhs->init (this->n_dofs(), this->n_local_dofs());
85 }
86 
87 
88 
89 void ExplicitSystem::assemble_qoi (const QoISet& qoi_indices)
90 {
91  // The user quantity of interest assembly gets to expect to
92  // accumulate on initially zero values
93  for (unsigned int i=0; i != qoi.size(); ++i)
94  if (qoi_indices.has_index(i))
95  qoi[i] = 0;
96 
97  Parent::assemble_qoi (qoi_indices);
98 }
99 
100 
101 
103 {
104  // The user quantity of interest derivative assembly gets to expect
105  // to accumulate on initially zero vectors
106  for (unsigned int i=0; i != qoi.size(); ++i)
107  if (qoi_indices.has_index(i))
108  this->add_adjoint_rhs(i).zero();
109 
110  Parent::assemble_qoi_derivative (qoi_indices);
111 }
112 
113 
114 
116 {
117  // Assemble the linear system
118  this->assemble ();
119 
120  // Update the system after the solve
121  this->update();
122 }
123 
124 
125 
127 {
128  // Possible that we cleared the _vectors but
129  // forgot to NULL-out the rhs?
130  if (this->n_vectors() == 0) rhs = NULL;
131 
132 
133  // Only need to add the rhs if it isn't there
134  // already!
135  if (rhs == NULL)
136  rhs = &(this->add_vector ("RHS Vector", false));
137 
139 }
140 
141 } // namespace libMesh

Site Created By: libMesh Developers
Last modified: February 07 2014 16:57:04 UTC

Hosted By:
SourceForge.net Logo