CppNoddy  0.92
Loading...
Searching...
No Matches
Functions
CppNoddy::Utility Namespace Reference

Some utility methods associated with CppNoddy containers. More...

Functions

DenseVector< double > uniform_node_vector (const double &lower, const double &upper, const std::size_t &N)
 Return a DENSE vector with the nodal points of a uniform mesh distributed between the upper/lower bounds as specified. More...
 
DenseVector< double > power_node_vector (const double &lower, const double &upper, const std::size_t &N, const double &power)
 Return a DENSE vector with the nodal points of a non-uniform mesh distributed between the upper/lower bounds as specified with more nodes clustered near lower or upper depending upon the differencee of the power from unity. More...
 
DenseVector< double > two_uniform_node_vector (const double &lower, const double &mid, const double &upper, const std::size_t &N1, const std::size_t &N2)
 Return a dense vector with two uniform distributions in two separate regions. More...
 
DenseVector< double > three_uniform_node_vector (const double &lower, const double &mid1, const double &mid2, const double &upper, const std::size_t &N1, const std::size_t &N2, const std::size_t &N3)
 Return a dense vector with two uniform distributions in two separate regions. More...
 
DenseVector< double > mid_weighted_node_vector (const double &lower, const double &upper, const std::size_t &N, const double &power)
 Return a dense vector of nodal positions with more nodes concentrated at the mid point of the range. More...
 
template<typename _Type >
void vels_from_streamfn_Cartesian (const TwoD_Node_Mesh< _Type > &source, TwoD_Node_Mesh< _Type > &uv)
 
template<typename _Type >
void vels_from_streamfn_Cartesian (const TwoD_Mapped_Node_Mesh< _Type > &source, TwoD_Mapped_Node_Mesh< _Type > &uv)
 
DenseMatrix< double > multiply (DenseMatrix< double > &A, DenseMatrix< double > &B)
 BLAS wrapper to do DOUBLE DENSE A_{MxK} * B_{KxN} = C_{MxN} Since this is a Fortran library, it assumes a column_major format, but CppNoddy uses row_major. More...
 
template<typename _Type >
_Type dot (const DenseVector< _Type > &X, const DenseVector< _Type > &Y)
 Templated dot product. More...
 
template<typename _Type >
int sgn (const _Type &a)
 
DenseVector< double > real (const DenseVector< D_complex > &X)
 Return a double DENSE vector containing the real part of a complex DENSE vector. More...
 
DenseVector< double > imag (const DenseVector< D_complex > &X)
 Return a double DENSE vector containing the imaginary part of a complex DENSE vector. More...
 
std::string stringify (const int &val)
 Return an integer value as a string - useful for file naming. More...
 
std::string stringify (const double &val, int p)
 Return a double value as a string - useful for file naming. More...
 
double max_abs_location (OneD_Node_Mesh< double > &mesh, unsigned var)
 
double max_abs_location_range (OneD_Node_Mesh< double > &mesh, unsigned var, double left, double right)
 

Detailed Description

Some utility methods associated with CppNoddy containers.

Function Documentation

◆ dot()

template<typename _Type >
_Type CppNoddy::Utility::dot ( const DenseVector< _Type > &  X,
const DenseVector< _Type > &  Y 
)

Templated dot product.

Parameters
XFirst dense vector
YSecond dense vector
Returns
The dot product

Definition at line 314 of file Utility.h.

314 {
315 if(X.size() != Y.size()) {
316 std::string problem;
317 problem = "The Utilities::dot method has been called \n";
318 problem += "with two unequal length vectors.";
319 throw ExceptionGeom(problem, X.size(), Y.size());
320 }
321 return inner_product(X.begin(), X.end(), Y.begin(), _Type(0.0));
322 }
CppNoddy::DenseVector::size
std::size_t size() const
A pass-thru definition to get the size of the vector.
Definition: DenseVector.h:330
CppNoddy::DenseVector::begin
elt_iter begin()
Pass through to the storage container.
Definition: DenseVector.h:82
CppNoddy::DenseVector::end
elt_iter end()
Pass through to the storage container.
Definition: DenseVector.h:102
CppNoddy::ExceptionGeom
An exception class to be thrown when a container of incorrect geometry used in any class/method.
Definition: Exceptions.h:47

References CppNoddy::DenseVector< _Type >::begin(), CppNoddy::DenseVector< _Type >::end(), and CppNoddy::DenseVector< _Type >::size().

Referenced by CppNoddy::ODE_BVP< _Type, _Xtype >::arclength_solve(), CppNoddy::Newton< _Type >::arclength_solve(), CppNoddy::PDE_IBVP< _Type >::assemble_matrix_problem(), CppNoddy::Equation_1matrix< _Type, double >::get_jacobian_of_matrix0_mult_vector(), and main().

◆ imag()

DenseVector< double > CppNoddy::Utility::imag ( const DenseVector< D_complex > &  X)

Return a double DENSE vector containing the imaginary part of a complex DENSE vector.

Parameters
XThe complex vector to take the imaginary part of

Definition at line 185 of file Utility.cpp.

185 {
186 DenseVector<double> temp(X.size(), 0.0);
187 for(std::size_t i = 0; i < X.size(); ++i) {
188 temp[ i ] = X[ i ].imag();
189 }
190 return temp;
191 }
CppNoddy::DenseVector
An DenseVector class – a dense vector object.
Definition: DenseVector.h:34

References CppNoddy::DenseVector< _Type >::size().

◆ max_abs_location()

double CppNoddy::Utility::max_abs_location ( OneD_Node_Mesh< double > &  mesh,
unsigned  var 
)

Definition at line 54 of file Utility.cpp.

54 {
55 double max(0.0);
56 std::size_t maxIndex(0);
57 // step through the nodes
58 for ( std::size_t node = 0; node < mesh.get_nnodes(); ++node ) {
59 //std::cout << "m_X[node]=" << m_X[node] << " left=" << left << " right=" << right << "\n";
60 if ( std::abs(mesh(node,var)) > max ) {
61 maxIndex = node;
62 max = std::abs( mesh(node,var) );
63 }
64 }
65 if ( ( maxIndex == 0 ) || ( maxIndex == mesh.get_nnodes()-1 ) ) {
66 std::cout << "[WARNING] MaxAbsLocationRange: maximumum absolute nodal value is first/last node. \n";
67 return mesh.coord( maxIndex );
68 }
69 double f1,f2,f3;
70 double x1,x2,x3;
71 f1 = std::abs(mesh(maxIndex-1,var));
72 f2 = std::abs(mesh(maxIndex,var));
73 f3 = std::abs(mesh(maxIndex+1,var));
74 x1 = mesh.coord(maxIndex-1);
75 x2 = mesh.coord(maxIndex);
76 x3 = mesh.coord(maxIndex+1);
77 return ( f1*(x2+x3)/((x1-x2)*(x1-x3)) + f2*(x1+x3)/((x2-x1)*(x2-x3)) + f3*(x1+x2)/((x3-x1)*(x3-x2)) )
78 / ( 2.*f1/((x1-x2)*(x1-x3)) + 2.*f2/((x2-x1)*(x2-x3)) + 2.*f3/((x3-x1)*(x3-x2)) );
79 }
CppNoddy::OneD_Node_Mesh::coord
const _Xtype & coord(const std::size_t &node) const
Access a nodal position.
Definition: OneD_Node_Mesh.h:278
CppNoddy::OneD_Node_Mesh::get_nnodes
std::size_t get_nnodes() const
Definition: OneD_Node_Mesh.cpp:50

References CppNoddy::OneD_Node_Mesh< _Type, _Xtype >::coord(), and CppNoddy::OneD_Node_Mesh< _Type, _Xtype >::get_nnodes().

◆ max_abs_location_range()

double CppNoddy::Utility::max_abs_location_range ( OneD_Node_Mesh< double > &  mesh,
unsigned  var,
double  left,
double  right 
)

Definition at line 82 of file Utility.cpp.

82 {
83 double max(0.0);
84 std::size_t maxIndex(0);
85 // step through the nodes
86 for ( std::size_t node = 0; node < mesh.get_nnodes(); ++node ) {
87 //std::cout << "m_X[node]=" << m_X[node] << " left=" << left << " right=" << right << "\n";
88 if ( ( mesh.coord(node) >= left ) && ( mesh.coord(node) <= right ) ) {
89 if ( std::abs(mesh(node,var)) > max ) {
90 maxIndex = node;
91 max = std::abs( mesh(node,var) );
92 }
93 }
94 }
95 if ( ( maxIndex == 0 ) || ( maxIndex == mesh.get_nnodes()-1 ) ) {
96 std::cout << "[WARNING] MaxAbsLocationRange: maximumum absolute nodal value is first/last node. \n";
97 return mesh.coord( maxIndex );
98 }
99 double f1,f2,f3;
100 double x1,x2,x3;
101 f1 = std::abs(mesh(maxIndex-1,var));
102 f2 = std::abs(mesh(maxIndex,var));
103 f3 = std::abs(mesh(maxIndex+1,var));
104 x1 = mesh.coord(maxIndex-1);
105 x2 = mesh.coord(maxIndex);
106 x3 = mesh.coord(maxIndex+1);
107 return ( f1*(x2+x3)/((x1-x2)*(x1-x3)) + f2*(x1+x3)/((x2-x1)*(x2-x3)) + f3*(x1+x2)/((x3-x1)*(x3-x2)) )
108 / ( 2.*f1/((x1-x2)*(x1-x3)) + 2.*f2/((x2-x1)*(x2-x3)) + 2.*f3/((x3-x1)*(x3-x2)) );
109 }

References CppNoddy::OneD_Node_Mesh< _Type, _Xtype >::coord(), and CppNoddy::OneD_Node_Mesh< _Type, _Xtype >::get_nnodes().

◆ mid_weighted_node_vector()

DenseVector< double > CppNoddy::Utility::mid_weighted_node_vector ( const double &  lower,
const double &  upper,
const std::size_t &  N,
const double &  power 
)

Return a dense vector of nodal positions with more nodes concentrated at the mid point of the range.

Parameters
lowerThe first nodal position.
upperThe final nodal position.
NThe number of nodes required.
powerA measure of the non-uniformity, power = 1 => uniform distribution

Definition at line 155 of file Utility.cpp.

155 {
156 DenseVector<double> node_vector(N, 0.0);
157 // make a center weighted distribution over -1 to 1
158 for(std::size_t i = 0; i < N; ++i) {
159 double s(-1.0 + 2.0 * i / (N - 1));
160 node_vector[ i ] = (weight * std::pow(s, 3) + s) / (weight + 1);
161 }
162 // map the -1 to 1 range to lower to upper
163 for(std::size_t i = 0; i < N; ++i) {
164 // move the range to 0->2
165 node_vector[ i ] += 1.0;
166 // move the range tp 0->1
167 node_vector[ i ] /= 2.0;
168 // move the range to lower -> upper
169 node_vector[ i ] *= (upper - lower);
170 node_vector[ i ] += lower;
171 }
172 return node_vector;
173 }

◆ multiply()

DenseMatrix< double > CppNoddy::Utility::multiply ( DenseMatrix< double > &  A,
DenseMatrix< double > &  B 
)

BLAS wrapper to do DOUBLE DENSE A_{MxK} * B_{KxN} = C_{MxN} Since this is a Fortran library, it assumes a column_major format, but CppNoddy uses row_major.

To be consistent we'll simply do (B^T)_{NxK} * (A^T)_{KxM} = (C^T)_{NxM} instead. Note that inversion of the transpose of the result C^T is handled implicitly via the construction of C.

Parameters
AFirst dense double matrix to be multiplied
BSecond dense double matrix to be multiplied
Returns
The result of the multiplication C=A*B

Definition at line 225 of file Utility.cpp.

225 {
226#ifndef LAPACK
227 std::string problem;
228 problem = "The Utilities::multiply method has been called\n";
229 problem += "but the compiler option -DLAPACK was not provided when\n";
230 problem += "the library was built. This non-member function requires BLAS.";
231 throw ExceptionExternal(problem);
232#else
233 // set the matrix geometries
234 std::size_t M = A.nrows();
235 std::size_t N = B.ncols();
236 std::size_t K = A.ncols();
237 // No need to transpose first, because we will in fact do B^T * A^T = C^T
238 // New the memory for the result
239 FortranData Af(A, false);
240 FortranData Bf(B, false);
241 FortranData Cf(M * N);
242#ifdef PARANOID
243
244 if(K != B.nrows()) {
245 std::string problem(" The LAPACK::multiply method has detected a failure \n");
246 throw ExceptionGeom(problem, A.nrows(), A.ncols(), B.nrows(), B.ncols());
247 }
248#endif
249 // call Fortran BLAS
250 // call Fortran BLAS
251 BLAS_DGEMM((char*) "N", (char*) "N", N, M, K, 1.0, Bf.base(), N, Af.base(), K, 0.0, Cf.base(), N);
252 // Return a DenseMatrix<double> from the results -- since Cf is in column_major
253 // format, this will actually transpose the Cf data to provide C as required
254 return (Cf.to_dense_matrix(M, N));
255#endif
256
257 }
BLAS_DGEMM
#define BLAS_DGEMM(TRANSA, TRANSB, M, N, K, ALPHA, A, LDA, B, LDB, BETA, C, LDC)
The BLAS DOUBLE dot-product Fortran interface : FUNCTION DDOT(N,X,INCX,Y,INCY)
Definition: FortranBLAS.h:19
CppNoddy::DenseMatrix::nrows
std::size_t nrows() const
Definition: DenseMatrix.h:375
CppNoddy::DenseMatrix::ncols
std::size_t ncols() const
Definition: DenseMatrix.h:380
CppNoddy::ExceptionExternal
An exception to indicate that an error has been detected in an external (LAPACK) routine.
Definition: Exceptions.h:20
CppNoddy::FortranData
A little (legacy) utility class for passing CppNoddy containers to Fortran library routines.
Definition: FortranData.h:13

References CppNoddy::FortranData::base(), BLAS_DGEMM, CppNoddy::DenseMatrix< _Type >::ncols(), CppNoddy::DenseMatrix< _Type >::nrows(), and CppNoddy::FortranData::to_dense_matrix().

Referenced by main().

◆ power_node_vector()

DenseVector< double > CppNoddy::Utility::power_node_vector ( const double &  lower,
const double &  upper,
const std::size_t &  N,
const double &  power 
)

Return a DENSE vector with the nodal points of a non-uniform mesh distributed between the upper/lower bounds as specified with more nodes clustered near lower or upper depending upon the differencee of the power from unity.

When power=1 this should provide a uniform mesh.

Parameters
lowerThe lower bound of the uniform nodal distribution
upperThe upper bound of the uniform nodal distribution
NThe number of nodal points
powerA measure of the non-uniformity
Returns
A vector of nodal positions with a power law distribution

Definition at line 123 of file Utility.cpp.

123 {
124 DenseVector<double> V(N, 0.0);
125 for(std::size_t i = 0; i < N; ++i) {
126 V[ i ] = lower + (upper - lower) * std::pow((double)i / (N - 1), power);
127 }
128 return V;
129 }
V
@ V
Definition: BVPKarman.cpp:20

References V.

Referenced by main().

◆ real()

DenseVector< double > CppNoddy::Utility::real ( const DenseVector< D_complex > &  X)

Return a double DENSE vector containing the real part of a complex DENSE vector.

Parameters
XThe complex vector to take the real part of

Definition at line 177 of file Utility.cpp.

177 {
178 DenseVector<double> temp(X.size(), 0.0);
179 for(std::size_t i = 0; i < X.size(); ++i) {
180 temp[ i ] = X[ i ].real();
181 }
182 return temp;
183 }

References CppNoddy::DenseVector< _Type >::size().

◆ sgn()

template<typename _Type >
int CppNoddy::Utility::sgn ( const _Type &  a)

Definition at line 325 of file Utility.h.

325 {
326 if(a > (_Type)0) {
327 return 1;
328 } else if(a < (_Type)0) {
329 return -1;
330 } else {
331 return 0;
332 }
333 }

◆ stringify() [1/2]

std::string CppNoddy::Utility::stringify ( const double &  val,
int  p 
)

Return a double value as a string - useful for file naming.

Parameters
valThe double value to be stringified
pPrecision to be used in the output

Definition at line 199 of file Utility.cpp.

199 {
200 std::stringstream temp;
201 temp.precision(p);
202 temp << val;
203 return temp.str();
204 }

References p.

◆ stringify() [2/2]

std::string CppNoddy::Utility::stringify ( const int &  val)

Return an integer value as a string - useful for file naming.

Parameters
valThe integer value to be stringified.

Definition at line 193 of file Utility.cpp.

193 {
194 std::stringstream temp;
195 temp << val;
196 return temp.str();
197 }

Referenced by CppNoddy::TrackerFile::header(), main(), and CppNoddy::TwoD_Node_Mesh< _Type >::remesh1().

◆ three_uniform_node_vector()

DenseVector< double > CppNoddy::Utility::three_uniform_node_vector ( const double &  lower,
const double &  mid1,
const double &  mid2,
const double &  upper,
const std::size_t &  N1,
const std::size_t &  N2,
const std::size_t &  N3 
)

Return a dense vector with two uniform distributions in two separate regions.

Parameters
lowerThe first node
mid1The node that defines the first interior boundary
mid2The node that defines the second interior boundary
upperThe final node
N1The number of nodes in the first region
N2The number of nodes in the second region
N3The number of nodes in the third region
Returns
A combined vector of nodes of length N1+N2+N3

Definition at line 141 of file Utility.cpp.

141 {
142 DenseVector<double> first = uniform_node_vector(lower, mid1, N1);
143 DenseVector<double> second = uniform_node_vector(mid1, mid2, N2+1);
144 DenseVector<double> third = uniform_node_vector(mid2, upper, N3+1);
145 // skip the common elt by starting at 1 not 0
146 for(std::size_t i = 1; i < N2+1; ++i) {
147 first.push_back(second[ i ]);
148 }
149 for(std::size_t i = 1; i < N3+1; ++i) {
150 first.push_back(third[ i ]);
151 }
152 return first;
153 }
CppNoddy::DenseVector::push_back
void push_back(const _Type &fill)
A pass-thru definition of push_back.
Definition: DenseVector.h:310
CppNoddy::Utility::uniform_node_vector
DenseVector< double > uniform_node_vector(const double &lower, const double &upper, const std::size_t &N)
Return a DENSE vector with the nodal points of a uniform mesh distributed between the upper/lower bou...
Definition: Utility.cpp:113

References CppNoddy::DenseVector< _Type >::push_back(), and uniform_node_vector().

◆ two_uniform_node_vector()

DenseVector< double > CppNoddy::Utility::two_uniform_node_vector ( const double &  lower,
const double &  mid,
const double &  upper,
const std::size_t &  N1,
const std::size_t &  N2 
)

Return a dense vector with two uniform distributions in two separate regions.

Parameters
lowerThe first node
midThe node that defines the boundary between the uniform meshes
upperThe final node
N1The number of nodes in the first region
N2The number of nodes in the second region
Returns
A combined vector of nodes of length N1+N2

Definition at line 131 of file Utility.cpp.

131 {
132 DenseVector<double> first = uniform_node_vector(lower, mid, N1);
133 DenseVector<double> second = uniform_node_vector(mid, upper, N2 + 1);
134 // skip the common elt by starting at 1 not 0
135 for(std::size_t i = 1; i < N2+1; ++i) {
136 first.push_back(second[ i ]);
137 }
138 return first;
139 }

References CppNoddy::DenseVector< _Type >::push_back(), and uniform_node_vector().

◆ uniform_node_vector()

DenseVector< double > CppNoddy::Utility::uniform_node_vector ( const double &  lower,
const double &  upper,
const std::size_t &  N 
)

Return a DENSE vector with the nodal points of a uniform mesh distributed between the upper/lower bounds as specified.

Parameters
lowerThe lower bound of the uniform nodal distribution
upperThe upper bound of the uniform nodal distribution
NThe number of nodal points

Definition at line 113 of file Utility.cpp.

113 {
114 DenseVector<double> V;
115 V.reserve(N);
116 const double delta = (upper - lower) / (N - 1);
117 for(std::size_t i = 0; i < N; ++i) {
118 V.push_back(lower + delta * i);
119 }
120 return V;
121 }

References V.

Referenced by CppNoddy::FnQuadrature::FnQuadrature(), CppNoddy::TwoD_Mapped_Node_Mesh< _Type >::init_mapping(), main(), CppNoddy::Example::Neutral_residual::Neutral_residual(), CppNoddy::FnQuadrature::set_subintervals(), three_uniform_node_vector(), and two_uniform_node_vector().

◆ vels_from_streamfn_Cartesian() [1/2]

template<typename _Type >
void CppNoddy::Utility::vels_from_streamfn_Cartesian ( const TwoD_Mapped_Node_Mesh< _Type > &  source,
TwoD_Mapped_Node_Mesh< _Type > &  uv 
)

Definition at line 166 of file Utility.h.

166 {
167 std::cout << "MAPPED MESH version\n";
168 std::size_t Nx(source.get_nnodes().first);
169 std::size_t Ny(source.get_nnodes().second);
170 double dX( source.get_comp_step_sizes().first );
171 double dY( source.get_comp_step_sizes().second );
172 // differentiate the streamfunction to get the velocity field
173 {
174 // west internal nodes
175 std::size_t i(0);
176 for(std::size_t j = 1; j < Ny - 1; ++j) {
177 double x = source.coord(i,j).first;
178 double Xd = source.FnComp_Xd(x);
179 double y = source.coord(i,j).second;
180 double Yd = source.FnComp_Yd(y);
181 uv(i, j, 0) = (source(i, j + 1, 0) - source(i, j - 1, 0)) * Yd / (2 * dY);
182 uv(i, j, 1) = -(-source(i + 2, j, 0) + 4. * source(i + 1, j, 0) - 3. * source(i, j, 0)) * Xd / (2 * dX);
183 }
184 }
185 {
186 // east internal nodes
187 std::size_t i(Nx - 1);
188 for(std::size_t j = 1; j < Ny - 1; ++j) {
189 double x = source.coord(i,j).first;
190 double Xd = source.FnComp_Xd(x);
191 double y = source.coord(i,j).second;
192 double Yd = source.FnComp_Yd(y);
193 uv(i, j, 0) = (source(i, j + 1, 0) - source(i, j - 1, 0)) * Yd / (2 * dY);
194 uv(i, j, 1) = -(source(i - 2, j, 0) - 4. * source(i - 1, j, 0) + 3. * source(i, j, 0)) * Xd / (2 * dX);
195 }
196 }
197 {
198 // south internal nodes
199 std::size_t j(0);
200 for(std::size_t i = 1; i < Nx - 1; ++i) {
201 double x = source.coord(i,j).first;
202 double Xd = source.FnComp_Xd(x);
203 double y = source.coord(i,j).second;
204 double Yd = source.FnComp_Yd(y);
205 uv(i, j, 0) = (-source(i, j + 2, 0) + 4. * source(i, j + 1, 0) - 3. * source(i, j, 0))* Yd / (2 * dY);
206 uv(i, j, 1) = -(source(i + 1, j, 0) - source(i - 1, j, 0)) * Xd / (2 * dX);
207 }
208 }
209 {
210 // north internal nodes
211 std::size_t j(Ny - 1);
212 for(std::size_t i = 1; i < Nx - 1; ++i) {
213 double x = source.coord(i,j).first;
214 double Xd = source.FnComp_Xd(x);
215 double y = source.coord(i,j).second;
216 double Yd = source.FnComp_Yd(y);
217 uv(i, j, 0) = (source(i, j - 2, 0) - 4. * source(i, j - 1, 0) + 3. * source(i, j, 0)) * Yd / (2 * dY);
218 uv(i, j, 1) = -(source(i + 1, j, 0) - source(i - 1, j, 0)) * Xd / (2 * dX);
219 }
220 }
221 {
222 // corner nodes
223 {
224 // sw
225 std::size_t i(0);
226 std::size_t j(0);
227 double x = source.coord(i,j).first;
228 double Xd = source.FnComp_Xd(x);
229 double y = source.coord(i,j).second;
230 double Yd = source.FnComp_Yd(y);
231 uv(i, j, 0) = (-source(i, j + 2, 0) + 4. * source(i, j + 1, 0) - 3. * source(i, j, 0))* Xd / (2 * dY);
232 uv(i, j, 1) = -(-source(i + 2, j, 0) + 4. * source(i + 1, j, 0) - 3. * source(i, j, 0))* Yd / (2 * dY);
233 }
234 {
235 // nw
236 std::size_t i(0);
237 std::size_t j(Ny - 1);
238 double x = source.coord(i,j).first;
239 double Xd = source.FnComp_Xd(x);
240 double y = source.coord(i,j).second;
241 double Yd = source.FnComp_Yd(y);
242 uv(i, j, 0) = (source(i, j - 2, 0) - 4. * source(i, j - 1, 0) + 3. * source(i, j, 0)) * Yd / (2 * dY);
243 uv(i, j, 1) = -(-source(i + 2, j, 0) + 4. * source(i + 1, j, 0) - 3. * source(i, j, 0))* Xd / (2 * dX);
244 }
245 {
246 // ne
247 std::size_t i(Nx - 1);
248 std::size_t j(Ny - 1);
249 double x = source.coord(i,j).first;
250 double Xd = source.FnComp_Xd(x);
251 double y = source.coord(i,j).second;
252 double Yd = source.FnComp_Yd(y);
253 uv(i, j, 0) = (source(i, j - 2, 0) - 4. * source(i, j - 1, 0) + 3. * source(i, j, 0)) * Yd / (2 * dY);
254 uv(i, j, 1) = -(source(i - 2, j, 0) - 4. * source(i - 1, j, 0) + 3. * source(i, j, 0)) * Xd/ (2 * dX);
255 }
256 {
257 // se
258 std::size_t i(Nx - 1);
259 std::size_t j(0);
260 double x = source.coord(i,j).first;
261 double Xd = source.FnComp_Xd(x);
262 double y = source.coord(i,j).second;
263 double Yd = source.FnComp_Yd(y);
264 uv(i, j, 0) = (-source(i, j + 2, 0) + 4. * source(i, j + 1, 0) - 3. * source(i, j, 0))* Yd / (2 * dY);
265 uv(i, j, 1) = -(source(i - 2, j, 0) - 4. * source(i - 1, j, 0) + 3. * source(i, j, 0)) * Xd / (2 * dX);
266 }
267 }
268 {
269 // interior nodes
270 for(std::size_t i = 1; i < Nx - 1; ++i) {
271 for(std::size_t j = 1; j < Ny - 1; ++j) {
272 double x = source.coord(i,j).first;
273 double Xd = source.FnComp_Xd(x);
274 double y = source.coord(i,j).second;
275 double Yd = source.FnComp_Yd(y);
276 uv(i, j, 0) = (source(i, j + 1, 0) - source(i, j - 1, 0)) * Yd / (2 * dY);
277 uv(i, j, 1) = -(source(i + 1, j, 0) - source(i - 1, j, 0)) * Xd / (2 * dX);
278 }
279 }
280 }
281 }
CppNoddy::Example::source
double source(const double &x, const double &y, const double &t)
Definition: IBVPLinear.cpp:26

◆ vels_from_streamfn_Cartesian() [2/2]

template<typename _Type >
void CppNoddy::Utility::vels_from_streamfn_Cartesian ( const TwoD_Node_Mesh< _Type > &  source,
TwoD_Node_Mesh< _Type > &  uv 
)

Definition at line 83 of file Utility.h.

83 {
84 std::cout << "PLAIN MESH version\n";
85 std::size_t Nx(source.get_nnodes().first);
86 std::size_t Ny(source.get_nnodes().second);
87 double dx(source.coord(1,1).first - source.coord(0,0).first);
88 double dy(source.coord(1,1).second - source.coord(0,0).second);
89 // differentiate the streamfunction to get the velocity field
90 {
91 // west internal nodes
92 std::size_t i(0);
93 for(std::size_t j = 1; j < Ny - 1; ++j) {
94 uv(i, j, 0) = (source(i, j + 1, 0) - source(i, j - 1, 0)) / (2 * dy);
95 uv(i, j, 1) = -(-source(i + 2, j, 0) + 4. * source(i + 1, j, 0) - 3. * source(i, j, 0)) / (2 * dx);
96 }
97 }
98 {
99 // east internal nodes
100 std::size_t i(Nx - 1);
101 for(std::size_t j = 1; j < Ny - 1; ++j) {
102 uv(i, j, 0) = (source(i, j + 1, 0) - source(i, j - 1, 0)) / (2 * dy);
103 uv(i, j, 1) = -(source(i - 2, j, 0) - 4. * source(i - 1, j, 0) + 3. * source(i, j, 0)) / (2 * dx);
104 }
105 }
106 {
107 // south internal nodes
108 std::size_t j(0);
109 for(std::size_t i = 1; i < Nx - 1; ++i) {
110 uv(i, j, 0) = (-source(i, j + 2, 0) + 4. * source(i, j + 1, 0) - 3. * source(i, j, 0)) / (2 * dy);
111 uv(i, j, 1) = -(source(i + 1, j, 0) - source(i - 1, j, 0)) / (2 * dx);
112 }
113 }
114 {
115 // north internal nodes
116 std::size_t j(Ny - 1);
117 for(std::size_t i = 1; i < Nx - 1; ++i) {
118 uv(i, j, 0) = (source(i, j - 2, 0) - 4. * source(i, j - 1, 0) + 3. * source(i, j, 0)) / (2 * dy);
119 uv(i, j, 1) = -(source(i + 1, j, 0) - source(i - 1, j, 0)) / (2 * dx);
120 }
121 }
122 {
123 // corner nodes
124 {
125 // sw
126 std::size_t i(0);
127 std::size_t j(0);
128 uv(i, j, 0) = (-source(i, j + 2, 0) + 4. * source(i, j + 1, 0) - 3. * source(i, j, 0)) / (2 * dy);
129 uv(i, j, 1) = -(-source(i + 2, j, 0) + 4. * source(i + 1, j, 0) - 3. * source(i, j, 0)) / (2 * dx);
130 }
131 {
132 // nw
133 std::size_t i(0);
134 std::size_t j(Ny - 1);
135 uv(i, j, 0) = (source(i, j - 2, 0) - 4. * source(i, j - 1, 0) + 3. * source(i, j, 0)) / (2 * dy);
136 uv(i, j, 1) = -(-source(i + 2, j, 0) + 4. * source(i + 1, j, 0) - 3. * source(i, j, 0)) / (2 * dx);
137 }
138 {
139 // ne
140 std::size_t i(Nx - 1);
141 std::size_t j(Ny - 1);
142 uv(i, j, 0) = (source(i, j - 2, 0) - 4. * source(i, j - 1, 0) + 3. * source(i, j, 0)) / (2 * dy);
143 uv(i, j, 1) = -(source(i - 2, j, 0) - 4. * source(i - 1, j, 0) + 3. * source(i, j, 0)) / (2 * dx);
144 }
145 {
146 // se
147 std::size_t i(Nx - 1);
148 std::size_t j(0);
149 uv(i, j, 0) = (-source(i, j + 2, 0) + 4. * source(i, j + 1, 0) - 3. * source(i, j, 0)) / (2 * dy);
150 uv(i, j, 1) = -(source(i - 2, j, 0) - 4. * source(i - 1, j, 0) + 3. * source(i, j, 0)) / (2 * dx);
151 }
152 }
153 {
154 // interior nodes
155 for(std::size_t i = 1; i < Nx - 1; ++i) {
156 for(std::size_t j = 1; j < Ny - 1; ++j) {
157 uv(i, j, 0) = (source(i, j + 1, 0) - source(i, j - 1, 0)) / (2 * dy);
158 uv(i, j, 1) = -(source(i + 1, j, 0) - source(i - 1, j, 0)) / (2 * dx);
159 }
160 }
161 }
162 }
Generated on Tue Jan 1 1980 00:00:00 for CppNoddy by   doxygen 1.9.5

© 2012

R.E. Hewitt