Nektar++
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
Public Member Functions | Static Public Member Functions | Static Public Attributes | Private Member Functions | Private Attributes | List of all members
Nektar::Utilities::NodeOpti2D3D Class Reference

#include <NodeOptiCAD.h>

Inheritance diagram for Nektar::Utilities::NodeOpti2D3D:
Inheritance graph
[legend]
Collaboration diagram for Nektar::Utilities::NodeOpti2D3D:
Collaboration graph
[legend]

Public Member Functions

 NodeOpti2D3D (NodeSharedPtr n, std::vector< ElUtilSharedPtr > e, ResidualSharedPtr r, std::map< LibUtilities::ShapeType, DerivUtilSharedPtr > d, optiType o, CADSurfSharedPtr s)
 
 ~NodeOpti2D3D ()
 
void Optimise ()
 
- Public Member Functions inherited from Nektar::Utilities::NodeOpti
 NodeOpti (NodeSharedPtr n, std::vector< ElUtilSharedPtr > e, ResidualSharedPtr r, std::map< LibUtilities::ShapeType, DerivUtilSharedPtr > d, optiType o, int dim)
 
virtual ~NodeOpti ()
 
void CalcMinJac ()
 
NodeOptiJobGetJob ()
 
template<int DIM>
NekDouble GetFunctional (NekDouble &minJacNew, bool gradient=true)
 Evaluate functional for elements connected to a node. More...
 
template<int DIM>
void MinEigen (NekDouble &val)
 Calculates minimum eigenvalue of Hessian matrix. More...
 
template<>
void MinEigen (NekDouble &val)
 
template<>
void MinEigen (NekDouble &val)
 

Static Public Member Functions

static NodeOptiSharedPtr create (NodeSharedPtr n, std::vector< ElUtilSharedPtr > e, ResidualSharedPtr r, std::map< LibUtilities::ShapeType, DerivUtilSharedPtr > d, optiType o)
 

Static Public Attributes

static int m_type
 

Private Member Functions

void ProcessGradient ()
 

Private Attributes

CADSurfSharedPtr surf
 

Additional Inherited Members

- Protected Member Functions inherited from Nektar::Utilities::NodeOpti
template<int DIM>
int IsIndefinite ()
 Returns 1 if Hessian matrix is indefinite and 0 otherwise. More...
 
template<>
int IsIndefinite ()
 
template<>
int IsIndefinite ()
 
- Static Protected Member Functions inherited from Nektar::Utilities::NodeOpti
static NekDouble c1 ()
 
static NekDouble gradTol ()
 
static NekDouble alphaTol ()
 
- Protected Attributes inherited from Nektar::Utilities::NodeOpti
NodeSharedPtr m_node
 
boost::mutex mtx
 
std::map
< LibUtilities::ShapeType,
std::vector< ElUtilSharedPtr > > 
m_data
 
Array< OneD, NekDoublem_grad
 
std::vector< NekDoublem_tmpStore
 
boost::unordered_map
< LibUtilities::ShapeType,
DerivArray
m_derivs
 
NekDouble m_minJac
 
ResidualSharedPtr m_res
 
std::map
< LibUtilities::ShapeType,
DerivUtilSharedPtr
m_derivUtils
 
optiType m_opti
 

Detailed Description

Definition at line 77 of file NodeOptiCAD.h.

Constructor & Destructor Documentation

Nektar::Utilities::NodeOpti2D3D::NodeOpti2D3D ( NodeSharedPtr  n,
std::vector< ElUtilSharedPtr e,
ResidualSharedPtr  r,
std::map< LibUtilities::ShapeType, DerivUtilSharedPtr d,
optiType  o,
CADSurfSharedPtr  s 
)
inline

Definition at line 80 of file NodeOptiCAD.h.

Referenced by create().

84  : NodeOpti(n, e, r, d, o, 3), surf(s)
85  {
86  }
NodeOpti(NodeSharedPtr n, std::vector< ElUtilSharedPtr > e, ResidualSharedPtr r, std::map< LibUtilities::ShapeType, DerivUtilSharedPtr > d, optiType o, int dim)
Definition: NodeOpti.h:59
Nektar::Utilities::NodeOpti2D3D::~NodeOpti2D3D ( )
inline

Definition at line 88 of file NodeOptiCAD.h.

88 {};

Member Function Documentation

static NodeOptiSharedPtr Nektar::Utilities::NodeOpti2D3D::create ( NodeSharedPtr  n,
std::vector< ElUtilSharedPtr e,
ResidualSharedPtr  r,
std::map< LibUtilities::ShapeType, DerivUtilSharedPtr d,
optiType  o 
)
inlinestatic

Definition at line 93 of file NodeOptiCAD.h.

References NodeOpti2D3D().

96  {
97  std::vector<std::pair<int, CADSurfSharedPtr> > ss = n->GetCADSurfs();
98  return NodeOptiSharedPtr(new NodeOpti2D3D(n, e, r, d, o, ss[0].second));
99  }
NodeOpti2D3D(NodeSharedPtr n, std::vector< ElUtilSharedPtr > e, ResidualSharedPtr r, std::map< LibUtilities::ShapeType, DerivUtilSharedPtr > d, optiType o, CADSurfSharedPtr s)
Definition: NodeOptiCAD.h:80
boost::shared_ptr< NodeOpti > NodeOptiSharedPtr
Definition: NodeOpti.h:135
void Nektar::Utilities::NodeOpti2D3D::Optimise ( )
virtual

Implements Nektar::Utilities::NodeOpti.

Definition at line 151 of file NodeOptiCAD.cpp.

References ASSERTL0, and CellMLToNektar.cellml_metadata::p.

152 {
153  NekDouble minJacNew;
154  NekDouble currentW = GetFunctional<3>(minJacNew);
155  NekDouble newVal = currentW;
156 
157  if (m_grad[0] * m_grad[0] + m_grad[1] * m_grad[1] + m_grad[2] * m_grad[2] >
158  gradTol())
159  {
160  // modify the gradient to be on the cad system
161  ProcessGradient();
162 
163  // needs to optimise
164  Array<OneD, NekDouble> uvc = m_node->GetCADSurfInfo(surf->GetId());
165  NekDouble xc = m_node->m_x;
166  NekDouble yc = m_node->m_y;
167  NekDouble zc = m_node->m_z;
168  Array<OneD, NekDouble> uvt(2);
169  Array<OneD, NekDouble> p;
170  Array<OneD, NekDouble> bd = surf->GetBounds();
171 
172  Array<OneD, NekDouble> sk(2);
173  NekDouble val;
174 
175  // Calculate minimum eigenvalue
176  MinEigen<2>(val);
177 
178  if (val < 1e-6)
179  {
180  // Add constant identity to Hessian matrix.
181  m_grad[2] += 1e-6 - val;
182  m_grad[4] += 1e-6 - val;
183  }
184 
185  sk[0] = -1.0 / (m_grad[2] * m_grad[4] - m_grad[3] * m_grad[3]) *
186  (m_grad[4] * m_grad[0] - m_grad[3] * m_grad[1]);
187  sk[1] = -1.0 / (m_grad[2] * m_grad[4] - m_grad[3] * m_grad[3]) *
188  (m_grad[2] * m_grad[1] - m_grad[3] * m_grad[0]);
189 
190  bool found = false;
191  NekDouble pg = (m_grad[0] * sk[0] + m_grad[1] * sk[1]);
192 
193  // normal gradient line Search
194  NekDouble alpha = 1.0;
195  while (alpha > alphaTol())
196  {
197  uvt[0] = uvc[0] + alpha * sk[0];
198  uvt[1] = uvc[1] + alpha * sk[1];
199 
200  if (uvt[0] < bd[0] || uvt[0] > bd[1] || uvt[1] < bd[2] ||
201  uvt[1] > bd[3])
202  {
203  alpha /= 2.0;
204  continue;
205  }
206 
207  p = surf->P(uvt);
208  m_node->m_x = p[0];
209  m_node->m_y = p[1];
210  m_node->m_z = p[2];
211 
212  newVal = GetFunctional<3>(minJacNew, false);
213 
214  // Wolfe conditions
215  if (newVal <= currentW + c1() * (alpha * pg))
216  {
217  found = true;
218  break;
219  }
220 
221  alpha /= 2.0;
222  }
223 
224  if (!found)
225  {
226  // reset the node
227  p = surf->P(uvc);
228  m_node->m_x = p[0];
229  m_node->m_y = p[1];
230  m_node->m_z = p[2];
231 
232  mtx.lock();
233  m_res->nReset[1]++;
234  mtx.unlock();
235  }
236  else
237  {
238  m_minJac = minJacNew;
239  m_node->Move(p, surf->GetId(), uvt);
240  }
241 
242  mtx.lock();
243  m_res->val = max(sqrt((m_node->m_x - xc) * (m_node->m_x - xc) +
244  (m_node->m_y - yc) * (m_node->m_y - yc) +
245  (m_node->m_z - zc) * (m_node->m_z - zc)),
246  m_res->val);
247  mtx.unlock();
248 
249  Array<OneD, NekDouble> uva = m_node->GetCADSurfInfo(surf->GetId());
250  if (uva[0] < bd[0] || uva[0] > bd[1] || uva[1] < bd[2] ||
251  uva[1] > bd[3])
252  {
253  ASSERTL0(false, "something when very wrong and the node finished "
254  "out of its parameter plane");
255  }
256  }
257 
258  mtx.lock();
259  m_res->func += newVal;
260  mtx.unlock();
261 }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:198
static NekDouble alphaTol()
Definition: NodeOpti.h:129
Array< OneD, NekDouble > m_grad
Definition: NodeOpti.h:109
static NekDouble c1()
Definition: NodeOpti.h:121
double NekDouble
ResidualSharedPtr m_res
Definition: NodeOpti.h:117
static NekDouble gradTol()
Definition: NodeOpti.h:125
void Nektar::Utilities::NodeOpti2D3D::ProcessGradient ( )
private

Definition at line 283 of file NodeOptiCAD.cpp.

284 {
285  Array<OneD, NekDouble> uvc = m_node->GetCADSurfInfo(surf->GetId());
286 
287  Array<OneD, NekDouble> grad = m_grad;
288  m_grad = Array<OneD, NekDouble>(5, 0.0);
289 
290  Array<OneD, NekDouble> d2 = surf->D2(uvc);
291  // r[0] x
292  // r[1] y
293  // r[2] z
294  // r[3] dx/du a
295  // r[4] dy/du b
296  // r[5] dz/du c
297  // r[6] dx/dv d
298  // r[7] dy/dv e
299  // r[8] dz/dv f
300  // r[9] d2x/du2
301  // r[10] d2y/du2
302  // r[11] d2z/du2
303  // r[12] d2x/dv2
304  // r[13] d2y/dv2
305  // r[14] d2z/dv2
306  // r[15] d2x/dudv
307  // r[16] d2y/dudv
308  // r[17] d2z/dudv
309  //
310  // grad[0] d/dx
311  // grad[1] d/dy
312  // grad[2] d/dz
313  // grad[3] d2/dx2
314  // grad[4] d2/dxdy
315  // grad[5] d2/dxdz
316  // grad[6] d2/dy2
317  // grad[7] d2/dydz
318  // grad[8] d2/dz2
319 
320  // Gradients
321  m_grad[0] = d2[3] * grad[0] + d2[4] * grad[1] + d2[5] * grad[2];
322  m_grad[1] = d2[6] * grad[0] + d2[7] * grad[1] + d2[8] * grad[2];
323 
324  m_grad[2] = grad[0] * d2[9] + grad[1] * d2[10] + grad[2] * d2[11] +
325  grad[3] * d2[3] * d2[3] + grad[6] * d2[4] * d2[4] +
326  grad[8] * d2[5] * d2[5] + 2.0 * grad[4] * d2[3] * d2[4] +
327  2.0 * grad[5] * d2[3] * d2[5] + 2.0 * grad[7] * d2[4] * d2[5];
328 
329  m_grad[4] = grad[0] * d2[12] + grad[1] * d2[13] + grad[2] * d2[14] +
330  grad[3] * d2[6] * d2[6] + grad[6] * d2[7] * d2[7] +
331  grad[8] * d2[8] * d2[8] + 2.0 * grad[4] * d2[6] * d2[7] +
332  2.0 * grad[5] * d2[6] * d2[8] + 2.0 * grad[7] * d2[7] * d2[8];
333 
334  m_grad[3] = grad[0] * d2[15] + grad[1] * d2[16] + grad[2] * d2[17] +
335  grad[3] * d2[3] * d2[6] + grad[6] * d2[4] * d2[7] +
336  grad[8] * d2[5] * d2[8] +
337  grad[4] * (d2[3] * d2[7] + d2[4] * d2[6]) +
338  grad[5] * (d2[3] * d2[8] + d2[5] * d2[6]) +
339  grad[7] * (d2[4] * d2[8] + d2[5] * d2[7]);
340 }
Array< OneD, NekDouble > m_grad
Definition: NodeOpti.h:109

Member Data Documentation

int Nektar::Utilities::NodeOpti2D3D::m_type
static
Initial value:

Definition at line 92 of file NodeOptiCAD.h.

CADSurfSharedPtr Nektar::Utilities::NodeOpti2D3D::surf
private

Definition at line 103 of file NodeOptiCAD.h.