fish_poisson_node_update.cc
Go to the documentation of this file.
1 //LIC// ====================================================================
2 //LIC// This file forms part of oomph-lib, the object-oriented,
3 //LIC// multi-physics finite-element library, available
4 //LIC// at http://www.oomph-lib.org.
5 //LIC//
6 //LIC// Version 1.0; svn revision $LastChangedRevision$
7 //LIC//
8 //LIC// $LastChangedDate$
9 //LIC//
10 //LIC// Copyright (C) 2006-2016 Matthias Heil and Andrew Hazel
11 //LIC//
12 //LIC// This library is free software; you can redistribute it and/or
13 //LIC// modify it under the terms of the GNU Lesser General Public
14 //LIC// License as published by the Free Software Foundation; either
15 //LIC// version 2.1 of the License, or (at your option) any later version.
16 //LIC//
17 //LIC// This library is distributed in the hope that it will be useful,
18 //LIC// but WITHOUT ANY WARRANTY; without even the implied warranty of
19 //LIC// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
20 //LIC// Lesser General Public License for more details.
21 //LIC//
22 //LIC// You should have received a copy of the GNU Lesser General Public
23 //LIC// License along with this library; if not, write to the Free Software
24 //LIC// Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA
25 //LIC// 02110-1301 USA.
26 //LIC//
27 //LIC// The authors may be contacted at oomph-lib@maths.man.ac.uk.
28 //LIC//
29 //LIC//====================================================================
30 // Driver for solution of 2D Poisson equation in fish-shaped domain with
31 // adaptivity and mesh updates
32 
33 // Generic oomph-lib headers
34 #include "generic.h"
35 
36 // The Poisson equations
37 #include "poisson.h"
38 
39 // The fish mesh
40 #include "meshes/fish_mesh.h"
41 
42 using namespace std;;
43 
44 using namespace oomph;
45 
46 //============ start_of_namespace=====================================
47 /// Namespace for const source term in Poisson equation
48 //====================================================================
49 namespace ConstSourceForPoisson
50 {
51 
52  /// Strength of source function: default value -1.0
53  double Strength=-1.0;
54 
55 /// Const source function
56  void get_source(const Vector<double>& x, double& source)
57  {
58  source = Strength;
59  }
60 
61 } // end of namespace
62 
63 
64 
65 
66 //======start_of_problem_class========================================
67 /// Refineable Poisson problem in fish-shaped domain.
68 /// Template parameter identifies the element type.
69 //====================================================================
70 template<class ELEMENT>
71 class RefineableFishPoissonProblem : public Problem
72 {
73 
74 public:
75 
76  /// Constructor
78 
79  /// Destructor: Empty
81 
82  /// Update the problem specs after solve (empty)
84 
85  /// Update the problem specs before solve (empty)
87 
88  /// \short Overloaded version of the problem's access function to
89  /// the mesh. Recasts the pointer to the base Mesh object to
90  /// the actual mesh type.
91  RefineableFishMesh<ELEMENT>* mesh_pt()
92  {
93  return dynamic_cast<RefineableFishMesh<ELEMENT>*>(Problem::mesh_pt());
94  }
95 
96  /// \short Doc the solution. Output directory and labels are specified
97  /// by DocInfo object
98  void doc_solution(DocInfo& doc_info);
99 
100 }; // end of problem class
101 
102 
103 
104 
105 
106 //===========start_of_constructor=========================================
107 /// Constructor for adaptive Poisson problem in fish-shaped
108 /// domain.
109 //========================================================================
110 template<class ELEMENT>
112 {
113 
114  // Build fish mesh -- this is a coarse base mesh consisting
115  // of four elements. We'll refine/adapt the mesh later.
116  Problem::mesh_pt()=new RefineableFishMesh<ELEMENT>;
117 
118  // Create/set error estimator
119  mesh_pt()->spatial_error_estimator_pt()=new Z2ErrorEstimator;
120 
121  // Set the boundary conditions for this problem: All nodes are
122  // free by default -- just pin the ones that have Dirichlet conditions
123  // here. Since the boundary values are never changed, we set
124  // them here rather than in actions_before_newton_solve().
125  unsigned num_bound = mesh_pt()->nboundary();
126  for(unsigned ibound=0;ibound<num_bound;ibound++)
127  {
128  unsigned num_nod= mesh_pt()->nboundary_node(ibound);
129  for (unsigned inod=0;inod<num_nod;inod++)
130  {
131  // Pin the single scalar value at this node
132  mesh_pt()->boundary_node_pt(ibound,inod)->pin(0);
133 
134  // Assign the homogenous boundary condition to the one
135  // and only nodal value
136  mesh_pt()->boundary_node_pt(ibound,inod)->set_value(0,0.0);
137  }
138  }
139 
140  // Loop over elements and set pointers to source function
141  unsigned n_element = mesh_pt()->nelement();
142  for(unsigned i=0;i<n_element;i++)
143  {
144  // Upcast from FiniteElement to the present element
145  ELEMENT *el_pt = dynamic_cast<ELEMENT*>(mesh_pt()->element_pt(i));
146 
147  //Set the source function pointer
148  el_pt->source_fct_pt() = &ConstSourceForPoisson::get_source;
149  }
150 
151  // Setup the equation numbering scheme
152  cout <<"Number of equations: " << assign_eqn_numbers() << std::endl;
153 
154 } // end of constructor
155 
156 
157 
158 
159 //=======start_of_doc=====================================================
160 /// Doc the solution in tecplot format.
161 //========================================================================
162 template<class ELEMENT>
164 {
165 
166  ofstream some_file;
167  char filename[100];
168 
169  // Number of plot points in each coordinate direction.
170  unsigned npts;
171  npts=5;
172 
173  // Output solution
174  sprintf(filename,"%s/soln%i.dat",doc_info.directory().c_str(),
175  doc_info.number());
176  some_file.open(filename);
177  mesh_pt()->output(some_file,npts);
178  some_file.close();
179 
180  // Output boundaries
181  sprintf(filename,"%s/boundaries%i.dat",doc_info.directory().c_str(),
182  doc_info.number());
183  some_file.open(filename);
184  mesh_pt()->output_boundaries(some_file);
185  some_file.close();
186 
187 } // end of doc
188 
189 
190 
191 
192 
193 
194 
195 //=================start_of_main==========================================
196 /// Demonstrate how to solve 2D Poisson problem in
197 /// fish-shaped domain with black-box mesh adaptation
198 /// and domain updates in response to changes in the domain
199 /// shape.
200 //========================================================================
201 int main()
202 {
203 
204  //Set up the problem with 9 node refineable Poisson elements
206 
207  // Setup labels for output
208  //------------------------
209  DocInfo doc_info;
210 
211  // Set output directory
212  doc_info.set_directory("RESLT");
213 
214  // Adjust the domain shape by changing the width of the fish
215  //----------------------------------------------------------
216  unsigned nstep=3;
217  for (unsigned i=0;i<nstep;i++)
218  {
219  // Get pointer to GeomObject that defines the position of the
220  // fish's back:
221  GeomObject* fish_back_pt=problem.mesh_pt()->fish_back_pt();
222 
223  // Recast to pointer to Circle object to get access to the member function
224  // that sets the y-position of the Circle's centre and decrease its
225  // value, making the fish narrower
226  dynamic_cast<Circle*>(fish_back_pt)->y_c()-=0.1;
227 
228  // Update the domain shape in response to the changes in its
229  // boundary
230  problem.mesh_pt()->node_update();
231 
232  // Solve the problem, allowing for up to two levels of refinement
233  problem.newton_solve(2);
234 
235  //Output solution
236  problem.doc_solution(doc_info);
237 
238  //Increment counter for solutions
239  doc_info.number()++;
240  }
241 
242 } // end of main
243 
void actions_after_newton_solve()
Update the problem specs after solve (empty)
double Strength
Strength of source function: default value -1.0.
virtual ~RefineableFishPoissonProblem()
Destructor: Empty.
void get_source(const Vector< double > &x, double &source)
Const source function.
void actions_before_newton_solve()
Update the problem specs before solve (empty)
void doc_solution(DocInfo &doc_info)
Doc the solution. Output directory and labels are specified by DocInfo object.
Namespace for const source term in Poisson equation.
RefineableFishMesh< ELEMENT > * mesh_pt()
Overloaded version of the problem&#39;s access function to the mesh. Recasts the pointer to the base Mesh...