geompack_mesh.template.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// Copyright (C) 2006-2022 Matthias Heil and Andrew Hazel
7// LIC//
8// LIC// This library is free software; you can redistribute it and/or
9// LIC// modify it under the terms of the GNU Lesser General Public
10// LIC// License as published by the Free Software Foundation; either
11// LIC// version 2.1 of the License, or (at your option) any later version.
12// LIC//
13// LIC// This library is distributed in the hope that it will be useful,
14// LIC// but WITHOUT ANY WARRANTY; without even the implied warranty of
15// LIC// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
16// LIC// Lesser General Public License for more details.
17// LIC//
18// LIC// You should have received a copy of the GNU Lesser General Public
19// LIC// License along with this library; if not, write to the Free Software
20// LIC// Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA
21// LIC// 02110-1301 USA.
22// LIC//
23// LIC// The authors may be contacted at oomph-lib@maths.man.ac.uk.
24// LIC//
25// LIC//====================================================================
26#ifndef OOMPH_GEOMPACK_MESH_TEMPLATE_CC
27#define OOMPH_GEOMPACK_MESH_TEMPLATE_CC
28
30
31
32namespace oomph
33{
34 //========================================================================
35 /// Quadrilateral mesh generator; Uses input from Geompack++.
36 /// See: http://members.shaw.ca/bjoe/
37 /// Currently only for four-noded quads -- extension to higher-order
38 /// quads should be trivial (see the corresponding classes for
39 /// triangular meshes).
40 //========================================================================
41 template<class ELEMENT>
43 TimeStepper* time_stepper_pt)
44 {
45 // Mesh can only be built with four-noded 2D Qelements.
46 MeshChecker::assert_geometric_element<QElementGeometricBase, ELEMENT>(2, 2);
47
48 // Create space for elements
49 unsigned nelem = Tmp_mesh_pt->nelement();
50 Element_pt.resize(nelem);
51
52 // Create space for nodes
53 unsigned nnod = Tmp_mesh_pt->nnode();
54 Node_pt.resize(nnod);
55
56 // Set number of boundaries
57 unsigned nbound = Tmp_mesh_pt->nboundary();
58 set_nboundary(nbound);
59
60 // Loop over elements in scaffold mesh, visit their nodes
61 for (unsigned e = 0; e < nelem; e++)
62 {
63 Element_pt[e] = new ELEMENT;
64 }
65
66
67 // At the moment we can only do four-noded quads
68 if (finite_element_pt(0)->nnode() != 4)
69 {
70 std::string error_message =
71 "GeompackQuadMesh can currently only deal with\n";
72 error_message += "four noded quads! \n";
73 error_message += "Generalisation to higher-order elements should be \n";
74 error_message += "straightforward but hasn't been done yet.\n";
75 error_message += "Do you want to volunteer? \n";
76
77 throw OomphLibError(
78 error_message, OOMPH_CURRENT_FUNCTION, OOMPH_EXCEPTION_LOCATION);
79 }
80
81
82 // In the first instance build all nodes from within all the elements
83 unsigned nnod_el = Tmp_mesh_pt->finite_element_pt(0)->nnode();
84
85 // Loop over elements in scaffold mesh, visit their nodes
86 for (unsigned e = 0; e < nelem; e++)
87 {
88 // Loop over all nodes in element
89 for (unsigned j = 0; j < nnod_el; j++)
90 {
91 // Create new node, using the NEW element's construct_node
92 // member function
93 finite_element_pt(e)->construct_node(j, time_stepper_pt);
94 }
95 }
96
97
98 std::map<Node*, unsigned> global_number;
99 unsigned global_count = 0;
100 // Loop over elements in scaffold mesh, visit their nodes
101 for (unsigned e = 0; e < nelem; e++)
102 {
103 // Loop over all nodes in element
104 for (unsigned j = 0; j < nnod_el; j++)
105 {
106 // Pointer to node in the scaffold mesh
107 Node* scaffold_node_pt = Tmp_mesh_pt->finite_element_pt(e)->node_pt(j);
108
109 // Get the (pseudo-)global node number in scaffold mesh
110 // (It's zero [=default] if not visited this one yet)
111 unsigned j_global = global_number[scaffold_node_pt];
112
113 // Haven't done this one yet
114 if (j_global == 0)
115 {
116 // Give it a number (not necessarily the global node
117 // number in the scaffold mesh -- we just need something
118 // to keep track...)
119 global_count++;
120 global_number[scaffold_node_pt] = global_count;
121
122 // Copy new node, created using the NEW element's construct_node
123 // function into global storage, using the same global
124 // node number that we've just associated with the
125 // corresponding node in the scaffold mesh
126 Node_pt[global_count - 1] = finite_element_pt(e)->node_pt(j);
127
128 // Assign coordinates
129 Node_pt[global_count - 1]->x(0) = scaffold_node_pt->x(0);
130 Node_pt[global_count - 1]->x(1) = scaffold_node_pt->x(1);
131
132
133 // Get pointer to set of mesh boundaries that this
134 // scaffold node occupies; NULL if the node is not on any boundary
135 std::set<unsigned>* boundaries_pt;
136 scaffold_node_pt->get_boundaries_pt(boundaries_pt);
137
138 // Loop over the mesh boundaries that the node in the scaffold mesh
139 // occupies and assign new node to the same ones.
140 if (boundaries_pt != 0)
141 {
142 // Convert it to a boundary node
143 this->convert_to_boundary_node(Node_pt[global_count - 1]);
144 // Add the node to the boundaries
145 for (std::set<unsigned>::iterator it = boundaries_pt->begin();
146 it != boundaries_pt->end();
147 ++it)
148 {
149 add_boundary_node(*it, Node_pt[global_count - 1]);
150 }
151 }
152 }
153 // This one has already been done: Kill it
154 else
155 {
156 // Kill it
157 delete finite_element_pt(e)->node_pt(j);
158
159 // Overwrite the element's pointer to local node by
160 // pointer to the already existing node -- identified by
161 // the number kept in the map
162 finite_element_pt(e)->node_pt(j) = Node_pt[j_global - 1];
163 }
164 }
165 }
166 }
167
168} // namespace oomph
169#endif
e
Definition: cfortran.h:571
void build_from_scaffold(TimeStepper *time_stepper_pt)
Build mesh from scaffold.
Nodes are derived from Data, but, in addition, have a definite (Eulerian) position in a space of a gi...
Definition: nodes.h:906
virtual void get_boundaries_pt(std::set< unsigned > *&boundaries_pt)
Return a pointer to set of mesh boundaries that this node occupies; this will be overloaded by Bounda...
Definition: nodes.h:1365
double & x(const unsigned &i)
Return the i-th nodal coordinate.
Definition: nodes.h:1060
An OomphLibError object which should be thrown when an run-time error is encountered....
////////////////////////////////////////////////////////////////////// //////////////////////////////...
Definition: timesteppers.h:231
std::string string(const unsigned &i)
Return the i-th string or "" if the relevant string hasn't been defined.
//////////////////////////////////////////////////////////////////// ////////////////////////////////...