MFEM  v4.0
Finite element discretization library
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
nurbs.hpp
Go to the documentation of this file.
1 // Copyright (c) 2010, Lawrence Livermore National Security, LLC. Produced at
2 // the Lawrence Livermore National Laboratory. LLNL-CODE-443211. All Rights
3 // reserved. See file COPYRIGHT for details.
4 //
5 // This file is part of the MFEM library. For more information and source code
6 // availability see http://mfem.org.
7 //
8 // MFEM is free software; you can redistribute it and/or modify it under the
9 // terms of the GNU Lesser General Public License (as published by the Free
10 // Software Foundation) version 2.1 dated February 1999.
11 
12 #ifndef MFEM_NURBS
13 #define MFEM_NURBS
14 
15 #include "../config/config.hpp"
16 #include "../general/table.hpp"
17 #include "../linalg/vector.hpp"
18 #include "element.hpp"
19 #include "mesh.hpp"
20 #ifdef MFEM_USE_MPI
21 #include "../general/communication.hpp"
22 #endif
23 #include <iostream>
24 
25 namespace mfem
26 {
27 
28 class GridFunction;
29 
30 
32 {
33 protected:
34  static const int MaxOrder;
35 
38 
39 public:
40  /// Create KnotVector
41  KnotVector() { }
42  KnotVector(std::istream &input);
43  KnotVector(int Order_, int NCP);
44  KnotVector(const KnotVector &kv) { (*this) = kv; }
45 
46  KnotVector &operator=(const KnotVector &kv);
47 
48  int GetNE() const { return NumOfElements; }
49  int GetNKS() const { return NumOfControlPoints - Order; }
50  int GetNCP() const { return NumOfControlPoints; }
51  int GetOrder() const { return Order; }
52  int Size() const { return knot.Size(); }
53 
54  /// Count the number of elements
55  void GetElements();
56 
57  bool isElement(int i) const { return (knot(Order+i) != knot(Order+i+1)); }
58 
59  double getKnotLocation(double xi, int ni) const
60  { return (xi*knot(ni+1) + (1. - xi)*knot(ni)); }
61 
62  int findKnotSpan(double u) const;
63 
64  void CalcShape (Vector &shape, int i, double xi) const;
65  void CalcDShape(Vector &grad, int i, double xi) const;
66 
67  void Difference(const KnotVector &kv, Vector &diff) const;
68  void UniformRefinement(Vector &newknots) const;
69  /** Return a new KnotVector with elevated degree by repeating the endpoints
70  of the knot vector. */
71  KnotVector *DegreeElevate(int t) const;
72 
73  void Flip();
74 
75  void Print(std::ostream &out) const;
76 
77  /// Destroys KnotVector
79 
80  double &operator[](int i) { return knot(i); }
81  const double &operator[](int i) const { return knot(i); }
82 };
83 
84 
86 {
87 protected:
88  int ni, nj, nk, Dim;
89  double *data; // the layout of data is: (Dim x ni x nj x nk)
90 
92 
93  int sd, nd;
94 
95  void swap(NURBSPatch *np);
96 
97  // Special B-NET access functions
98  int SetLoopDirection(int dir);
99  inline double &operator()(int i, int j);
100  inline const double &operator()(int i, int j) const;
101 
102  void init(int dim_);
103 
104  NURBSPatch(NURBSPatch *parent, int dir, int Order, int NCP);
105 
106 public:
107  NURBSPatch(const NURBSPatch &orig);
108  NURBSPatch(std::istream &input);
109  NURBSPatch(const KnotVector *kv0, const KnotVector *kv1, int dim_);
110  NURBSPatch(const KnotVector *kv0, const KnotVector *kv1,
111  const KnotVector *kv2, int dim_);
113 
114  ~NURBSPatch();
115 
116  void Print(std::ostream &out) const;
117 
118  void DegreeElevate(int dir, int t);
119  void KnotInsert (int dir, const KnotVector &knot);
120  void KnotInsert (int dir, const Vector &knot);
121 
122  void KnotInsert(Array<KnotVector *> &knot);
123  void DegreeElevate(int t);
124  void UniformRefinement();
125 
126  // Return the number of components stored in the NURBSPatch
127  int GetNC() const { return Dim; }
128  int GetNKV() const { return kv.Size(); }
129  KnotVector *GetKV(int i) { return kv[i]; }
130 
131  // Standard B-NET access functions
132  inline double &operator()(int i, int j, int l);
133  inline const double &operator()(int i, int j, int l) const;
134 
135  inline double &operator()(int i, int j, int k, int l);
136  inline const double &operator()(int i, int j, int k, int l) const;
137 
138  static void Get3DRotationMatrix(double n[], double angle, double r,
139  DenseMatrix &T);
140  void FlipDirection(int dir);
141  void SwapDirections(int dir1, int dir2);
142  void Rotate3D(double normal[], double angle);
143  int MakeUniformDegree(int degree = -1);
144  friend NURBSPatch *Interpolate(NURBSPatch &p1, NURBSPatch &p2);
145  friend NURBSPatch *Revolve3D(NURBSPatch &patch, double n[], double ang,
146  int times);
147 };
148 
149 
150 #ifdef MFEM_USE_MPI
151 class ParNURBSExtension;
152 #endif
153 
154 class NURBSPatchMap;
155 
156 
158 {
159 #ifdef MFEM_USE_MPI
160  friend class ParNURBSExtension;
161 #endif
162  friend class NURBSPatchMap;
163 
164 protected:
165  int mOrder; // see GetOrder() for description
168  // global entity counts
170  // local entity counts
173 
174  Array<int> activeVert; // activeVert[glob_vert] = loc_vert or -1
177  Array<int> activeDof; // activeDof[glob_dof] = loc_dof + 1 or 0
178 
180  int own_topo;
184 
185  // periodic BC info:
186  // - dof 2 dof map
187  // - master and slave boundary indices
191 
192  // global offsets, meshOffsets == meshVertexOffsets
197 
198  // global offsets, spaceOffsets == dofOffsets
203 
205 
208  Array2D<int> el_to_IJK; // IJK are "knot-span" indices!
209  Array2D<int> bel_to_IJK; // they are NOT element indices!
210 
212 
213  inline int KnotInd(int edge) const;
214  inline KnotVector *KnotVec(int edge);
215  inline const KnotVector *KnotVec(int edge) const;
216  inline const KnotVector *KnotVec(int edge, int oedge, int *okv) const;
217 
218  void CheckPatches();
219  void CheckBdrPatches();
220 
221  void GetPatchKnotVectors (int p, Array<KnotVector *> &kv);
222  void GetPatchKnotVectors (int p, Array<const KnotVector *> &kv) const;
224  void GetBdrPatchKnotVectors(int p, Array<const KnotVector *> &kv) const;
225 
226  void SetOrderFromOrders();
228 
229  // periodic BC helper functions
230  void InitDofMap();
231  void ConnectBoundaries();
232  void ConnectBoundaries2D(int bnd0, int bnd1);
233  void ConnectBoundaries3D(int bnd0, int bnd1);
234  int DofMap(int dof) const
235  {
236  return (d_to_d.Size() > 0 )? d_to_d[dof] : dof;
237  };
238 
239  // also count the global NumOfVertices and the global NumOfDofs
240  void GenerateOffsets();
241  // count the global NumOfElements
242  void CountElements();
243  // count the global NumOfBdrElements
244  void CountBdrElements();
245 
246  // generate the mesh elements
247  void Get2DElementTopo(Array<Element *> &elements) const;
248  void Get3DElementTopo(Array<Element *> &elements) const;
249 
250  // generate the boundary mesh elements
251  void Get2DBdrElementTopo(Array<Element *> &boundary) const;
252  void Get3DBdrElementTopo(Array<Element *> &boundary) const;
253 
254 
255  // FE space generation functions
256 
257  // based on activeElem, count NumOfActiveDofs, generate el_dof,
258  // el_to_patch, el_to_IJK, activeDof map (global-to-local)
260 
261  // generate elem_to_global-dof table for the active elements
262  // define el_to_patch, el_to_IJK, activeDof (as bool)
265 
266  // call after GenerateElementDofTable
268 
269  // generate the bdr-elem_to_global-dof table for the active bdr. elements
270  // define bel_to_patch, bel_to_IJK
273 
274  // FE --> Patch translation functions
275  void GetPatchNets (const Vector &Nodes, int vdim);
276  void Get2DPatchNets(const Vector &Nodes, int vdim);
277  void Get3DPatchNets(const Vector &Nodes, int vdim);
278 
279  // Patch --> FE translation functions
280  // Side effects: delete the patches, update the weights from the patches
281  void SetSolutionVector (Vector &Nodes, int vdim);
282  void Set2DSolutionVector(Vector &Nodes, int vdim);
283  void Set3DSolutionVector(Vector &Nodes, int vdim);
284 
285  // determine activeVert, NumOfActiveVertices from the activeElem array
286  void GenerateActiveVertices();
287 
288  // determine activeBdrElem, NumOfActiveBdrElems
289  void GenerateActiveBdrElems();
290 
291  void MergeWeights(Mesh *mesh_array[], int num_pieces);
292 
293  // to be used by ParNURBSExtension constructor(s)
295 
296 public:
297  /// Copy constructor: deep copy
298  NURBSExtension(const NURBSExtension &orig);
299  /// Read-in a NURBSExtension
300  NURBSExtension(std::istream &input);
301  /** @brief Create a NURBSExtension with elevated order by repeating the
302  endpoints of the knot vectors and using uniform weights of 1. */
303  /** If a knot vector in @a parent already has order greater than or equal to
304  @a newOrder, it will be used unmodified. */
305  NURBSExtension(NURBSExtension *parent, int newOrder);
306  /** @brief Create a NURBSExtension with elevated knot vector orders (by
307  repeating the endpoints of the knot vectors and using uniform weights of
308  1) as given by the array @a newOrders. */
309  /** If a knot vector in @a parent already has order greater than or equal to
310  the corresponding entry in @a newOrder, it will be used unmodified. */
311  NURBSExtension(NURBSExtension *parent, const Array<int> &newOrders);
312  /// Construct a NURBSExtension by merging a partitioned NURBS mesh
313  NURBSExtension(Mesh *mesh_array[], int num_pieces);
314 
315  // Generate connections between boundaries, such as periodic BCs
317  const Array<int> &GetMaster() const { return master; };
318  Array<int> &GetMaster() { return master; };
319  const Array<int> &GetSlave() const { return slave; };
320  Array<int> &GetSlave() { return slave; };
321  void MergeGridFunctions(GridFunction *gf_array[], int num_pieces,
322  GridFunction &merged);
323 
324  /// Destroy a NURBSExtension
325  virtual ~NURBSExtension();
326 
327  // Print functions
328  void Print(std::ostream &out) const;
329  void PrintCharacteristics(std::ostream &out) const;
330 
331  // Meta data functions
332  int Dimension() const { return patchTopo->Dimension(); }
333  int GetNP() const { return patchTopo->GetNE(); }
334  int GetNBP() const { return patchTopo->GetNBE(); }
335 
336  /// Read-only access to the orders of all knot vectors.
337  const Array<int> &GetOrders() const { return mOrders; }
338  /** @brief If all orders are identical, return that number. Otherwise, return
339  NURBSFECollection::VariableOrder. */
340  int GetOrder() const { return mOrder; }
341 
342  int GetNKV() const { return NumOfKnotVectors; }
343 
344  int GetGNV() const { return NumOfVertices; }
345  int GetNV() const { return NumOfActiveVertices; }
346  int GetGNE() const { return NumOfElements; }
347  int GetNE() const { return NumOfActiveElems; }
348  int GetGNBE() const { return NumOfBdrElements; }
349  int GetNBE() const { return NumOfActiveBdrElems; }
350 
351  int GetNTotalDof() const { return NumOfDofs; }
352  int GetNDof() const { return NumOfActiveDofs; }
353 
354  // Knotvector read-only access function
355  const KnotVector *GetKnotVector(int i) const { return knotVectors[i]; }
356 
357  // Mesh generation functions
358  void GetElementTopo (Array<Element *> &elements) const;
359  void GetBdrElementTopo(Array<Element *> &boundary) const;
360 
361  bool HavePatches() const { return (patches.Size() != 0); }
362 
365 
366  void GetVertexLocalToGlobal(Array<int> &lvert_vert);
367  void GetElementLocalToGlobal(Array<int> &lelem_elem);
368 
369  // Load functions
370  void LoadFE(int i, const FiniteElement *FE) const;
371  void LoadBE(int i, const FiniteElement *BE) const;
372 
373  const Vector &GetWeights() const { return weights; }
374  Vector &GetWeights() { return weights; }
375 
376  // Translation functions: from FE coordinates to IJK patch
377  // format and vice versa
378  void ConvertToPatches(const Vector &Nodes);
379  void SetKnotsFromPatches();
380  void SetCoordsFromPatches(Vector &Nodes);
381 
382  // Read a GridFunction written patch-by-patch, e.g. with PrintSolution().
383  void LoadSolution(std::istream &input, GridFunction &sol) const;
384  // Write a GridFunction patch-by-patch.
385  void PrintSolution(const GridFunction &sol, std::ostream &out) const;
386 
387  // Refinement methods
388  // new_degree = max(old_degree, min(old_degree + rel_degree, degree))
389  void DegreeElevate(int rel_degree, int degree = 16);
390  void UniformRefinement();
392 };
393 
394 
395 #ifdef MFEM_USE_MPI
397 {
398 private:
399  int *partitioning;
400 
401  Table *GetGlobalElementDofTable();
402  Table *Get2DGlobalElementDofTable();
403  Table *Get3DGlobalElementDofTable();
404 
405  void SetActive(const int *partitioning, const Array<bool> &active_bel);
406  void BuildGroups(const int *partitioning, const Table &elem_dof);
407 
408 public:
410 
412 
414 
415  ParNURBSExtension(MPI_Comm comm, NURBSExtension *parent, int *partitioning,
416  const Array<bool> &active_bel);
417 
418  // Create a parallel version of 'parent' with partitioning as in
419  // 'par_parent'; the 'parent' object is destroyed.
420  // The 'parent' can be either a local NURBSExtension or a global one.
422  const ParNURBSExtension *par_parent);
423 
424  virtual ~ParNURBSExtension() { delete [] partitioning; }
425 };
426 #endif
427 
428 
430 {
431 private:
432  const NURBSExtension *Ext;
433 
434  int I, J, K, pOffset, opatch;
435  Array<int> verts, edges, faces, oedge, oface;
436 
437  inline static int F(const int n, const int N)
438  { return (n < 0) ? 0 : ((n >= N) ? 2 : 1); }
439 
440  inline static int Or1D(const int n, const int N, const int Or)
441  { return (Or > 0) ? n : (N - 1 - n); }
442 
443  inline static int Or2D(const int n1, const int n2,
444  const int N1, const int N2, const int Or);
445 
446  // also set verts, edges, faces, orientations etc
447  void GetPatchKnotVectors (int p, const KnotVector *kv[]);
448  void GetBdrPatchKnotVectors(int p, const KnotVector *kv[], int *okv);
449 
450 public:
451  NURBSPatchMap(const NURBSExtension *ext) { Ext = ext; }
452 
453  int nx() { return I + 1; }
454  int ny() { return J + 1; }
455  int nz() { return K + 1; }
456 
457  void SetPatchVertexMap(int p, const KnotVector *kv[]);
458  void SetPatchDofMap (int p, const KnotVector *kv[]);
459 
460  void SetBdrPatchVertexMap(int p, const KnotVector *kv[], int *okv);
461  void SetBdrPatchDofMap (int p, const KnotVector *kv[], int *okv);
462 
463  inline int operator()(const int i) const;
464  inline int operator[](const int i) const { return (*this)(i); }
465 
466  inline int operator()(const int i, const int j) const;
467 
468  inline int operator()(const int i, const int j, const int k) const;
469 };
470 
471 
472 // Inline function implementations
473 
474 inline double &NURBSPatch::operator()(int i, int j)
475 {
476  return data[j%sd + sd*(i + (j/sd)*nd)];
477 }
478 
479 inline const double &NURBSPatch::operator()(int i, int j) const
480 {
481  return data[j%sd + sd*(i + (j/sd)*nd)];
482 }
483 
484 inline double &NURBSPatch::operator()(int i, int j, int l)
485 {
486 #ifdef MFEM_DEBUG
487  if (data == 0 || i < 0 || i >= ni || j < 0 || j >= nj || nk > 0 ||
488  l < 0 || l >= Dim)
489  {
490  mfem_error("NURBSPatch::operator() 2D");
491  }
492 #endif
493 
494  return data[(i+j*ni)*Dim+l];
495 }
496 
497 inline const double &NURBSPatch::operator()(int i, int j, int l) const
498 {
499 #ifdef MFEM_DEBUG
500  if (data == 0 || i < 0 || i >= ni || j < 0 || j >= nj || nk > 0 ||
501  l < 0 || l >= Dim)
502  {
503  mfem_error("NURBSPatch::operator() const 2D");
504  }
505 #endif
506 
507  return data[(i+j*ni)*Dim+l];
508 }
509 
510 inline double &NURBSPatch::operator()(int i, int j, int k, int l)
511 {
512 #ifdef MFEM_DEBUG
513  if (data == 0 || i < 0 || i >= ni || j < 0 || j >= nj || k < 0 ||
514  k >= nk || l < 0 || l >= Dim)
515  {
516  mfem_error("NURBSPatch::operator() 3D");
517  }
518 #endif
519 
520  return data[(i+(j+k*nj)*ni)*Dim+l];
521 }
522 
523 inline const double &NURBSPatch::operator()(int i, int j, int k, int l) const
524 {
525 #ifdef MFEM_DEBUG
526  if (data == 0 || i < 0 || i >= ni || j < 0 || j >= nj || k < 0 ||
527  k >= nk || l < 0 || l >= Dim)
528  {
529  mfem_error("NURBSPatch::operator() const 3D");
530  }
531 #endif
532 
533  return data[(i+(j+k*nj)*ni)*Dim+l];
534 }
535 
536 
537 inline int NURBSExtension::KnotInd(int edge) const
538 {
539  int kv = edge_to_knot[edge];
540  return (kv >= 0) ? kv : (-1-kv);
541 }
542 
544 {
545  return knotVectors[KnotInd(edge)];
546 }
547 
548 inline const KnotVector *NURBSExtension::KnotVec(int edge) const
549 {
550  return knotVectors[KnotInd(edge)];
551 }
552 
553 inline const KnotVector *NURBSExtension::KnotVec(int edge, int oedge, int *okv)
554 const
555 {
556  int kv = edge_to_knot[edge];
557  if (kv >= 0)
558  {
559  *okv = oedge;
560  return knotVectors[kv];
561  }
562  else
563  {
564  *okv = -oedge;
565  return knotVectors[-1-kv];
566  }
567 }
568 
569 
570 // static method
571 inline int NURBSPatchMap::Or2D(const int n1, const int n2,
572  const int N1, const int N2, const int Or)
573 {
574  // Needs testing
575  switch (Or)
576  {
577  case 0: return n1 + n2*N1;
578  case 1: return n2 + n1*N2;
579  case 2: return n2 + (N1 - 1 - n1)*N2;
580  case 3: return (N1 - 1 - n1) + n2*N1;
581  case 4: return (N1 - 1 - n1) + (N2 - 1 - n2)*N1;
582  case 5: return (N2 - 1 - n2) + (N1 - 1 - n1)*N2;
583  case 6: return (N2 - 1 - n2) + n1*N2;
584  case 7: return n1 + (N2 - 1 - n2)*N1;
585  }
586 #ifdef MFEM_DEBUG
587  mfem_error("NURBSPatchMap::Or2D");
588 #endif
589  return -1;
590 }
591 
592 inline int NURBSPatchMap::operator()(const int i) const
593 {
594  const int i1 = i - 1;
595  switch (F(i1, I))
596  {
597  case 0: return verts[0];
598  case 1: return pOffset + Or1D(i1, I, opatch);
599  case 2: return verts[1];
600  }
601 #ifdef MFEM_DEBUG
602  mfem_error("NURBSPatchMap::operator() const 1D");
603 #endif
604  return -1;
605 }
606 
607 inline int NURBSPatchMap::operator()(const int i, const int j) const
608 {
609  const int i1 = i - 1, j1 = j - 1;
610  switch (3*F(j1, J) + F(i1, I))
611  {
612  case 0: return verts[0];
613  case 1: return edges[0] + Or1D(i1, I, oedge[0]);
614  case 2: return verts[1];
615  case 3: return edges[3] + Or1D(j1, J, -oedge[3]);
616  case 4: return pOffset + Or2D(i1, j1, I, J, opatch);
617  case 5: return edges[1] + Or1D(j1, J, oedge[1]);
618  case 6: return verts[3];
619  case 7: return edges[2] + Or1D(i1, I, -oedge[2]);
620  case 8: return verts[2];
621  }
622 #ifdef MFEM_DEBUG
623  mfem_error("NURBSPatchMap::operator() const 2D");
624 #endif
625  return -1;
626 }
627 
628 inline int NURBSPatchMap::operator()(const int i, const int j, const int k)
629 const
630 {
631  // Needs testing
632  const int i1 = i - 1, j1 = j - 1, k1 = k - 1;
633  switch (3*(3*F(k1, K) + F(j1, J)) + F(i1, I))
634  {
635  case 0: return verts[0];
636  case 1: return edges[0] + Or1D(i1, I, oedge[0]);
637  case 2: return verts[1];
638  case 3: return edges[3] + Or1D(j1, J, oedge[3]);
639  case 4: return faces[0] + Or2D(i1, J - 1 - j1, I, J, oface[0]);
640  case 5: return edges[1] + Or1D(j1, J, oedge[1]);
641  case 6: return verts[3];
642  case 7: return edges[2] + Or1D(i1, I, oedge[2]);
643  case 8: return verts[2];
644  case 9: return edges[8] + Or1D(k1, K, oedge[8]);
645  case 10: return faces[1] + Or2D(i1, k1, I, K, oface[1]);
646  case 11: return edges[9] + Or1D(k1, K, oedge[9]);
647  case 12: return faces[4] + Or2D(J - 1 - j1, k1, J, K, oface[4]);
648  case 13: return pOffset + I*(J*k1 + j1) + i1;
649  case 14: return faces[2] + Or2D(j1, k1, J, K, oface[2]);
650  case 15: return edges[11] + Or1D(k1, K, oedge[11]);
651  case 16: return faces[3] + Or2D(I - 1 - i1, k1, I, K, oface[3]);
652  case 17: return edges[10] + Or1D(k1, K, oedge[10]);
653  case 18: return verts[4];
654  case 19: return edges[4] + Or1D(i1, I, oedge[4]);
655  case 20: return verts[5];
656  case 21: return edges[7] + Or1D(j1, J, oedge[7]);
657  case 22: return faces[5] + Or2D(i1, j1, I, J, oface[5]);
658  case 23: return edges[5] + Or1D(j1, J, oedge[5]);
659  case 24: return verts[7];
660  case 25: return edges[6] + Or1D(i1, I, oedge[6]);
661  case 26: return verts[6];
662  }
663 #ifdef MFEM_DEBUG
664  mfem_error("NURBSPatchMap::operator() const 3D");
665 #endif
666  return -1;
667 }
668 
669 }
670 
671 #endif
Array< KnotVector * > knotVectors
Definition: nurbs.hpp:182
Abstract class for Finite Elements.
Definition: fe.hpp:229
Array< KnotVector * > kv
Definition: nurbs.hpp:91
Array< int > f_meshOffsets
Definition: nurbs.hpp:195
KnotVector * DegreeElevate(int t) const
Definition: nurbs.cpp:57
NURBSPatch(NURBSPatch *parent, int dir, int Order, int NCP)
Definition: nurbs.cpp:404
Array2D< int > bel_to_IJK
Definition: nurbs.hpp:209
int Size() const
Logical size of the array.
Definition: array.hpp:118
void SetCoordsFromPatches(Vector &Nodes)
Definition: nurbs.cpp:2811
int Size() const
Definition: nurbs.hpp:52
int GetGNE() const
Definition: nurbs.hpp:346
Array< int > activeVert
Definition: nurbs.hpp:174
void ConnectBoundaries2D(int bnd0, int bnd1)
Definition: nurbs.cpp:1669
KnotVector * KnotVec(int edge)
Definition: nurbs.hpp:543
const Vector & GetWeights() const
Definition: nurbs.hpp:373
void swap(NURBSPatch *np)
Definition: nurbs.cpp:419
int GetNV() const
Definition: nurbs.hpp:345
const Array< int > & GetMaster() const
Definition: nurbs.hpp:317
Class for grid function - Vector with associated FE space.
Definition: gridfunc.hpp:27
void DegreeElevate(int dir, int t)
Definition: nurbs.cpp:697
int NumOfElements
Definition: nurbs.hpp:37
void LoadBE(int i, const FiniteElement *BE) const
Definition: nurbs.cpp:2779
void Generate2DElementDofTable()
Definition: nurbs.cpp:2498
int KnotInd(int edge) const
Definition: nurbs.hpp:537
void MergeGridFunctions(GridFunction *gf_array[], int num_pieces, GridFunction &merged)
Definition: nurbs.cpp:1960
Array< int > ldof_group
Definition: nurbs.hpp:411
Array< int > master
Definition: nurbs.hpp:189
virtual ~ParNURBSExtension()
Definition: nurbs.hpp:424
void Generate3DElementDofTable()
Definition: nurbs.cpp:2552
void KnotInsert(int dir, const KnotVector &knot)
Definition: nurbs.cpp:559
void MergeWeights(Mesh *mesh_array[], int num_pieces)
Definition: nurbs.cpp:1935
int GetNBE() const
Returns number of boundary elements.
Definition: mesh.hpp:679
void Get3DBdrElementTopo(Array< Element * > &boundary) const
Definition: nurbs.cpp:2430
NURBSPatchMap(const NURBSExtension *ext)
Definition: nurbs.hpp:451
void UniformRefinement()
Definition: nurbs.cpp:2950
int GetNE() const
Definition: nurbs.hpp:48
void Rotate3D(double normal[], double angle)
Definition: nurbs.cpp:1024
Array< int > e_meshOffsets
Definition: nurbs.hpp:194
int GetNKV() const
Definition: nurbs.hpp:342
int GetNKS() const
Definition: nurbs.hpp:49
int operator()(const int i) const
Definition: nurbs.hpp:592
void GetElements()
Count the number of elements.
Definition: nurbs.cpp:100
void SetOrdersFromKnotVectors()
Definition: nurbs.cpp:2153
Data type dense matrix using column-major storage.
Definition: densemat.hpp:23
int Size() const
Returns the size of the vector.
Definition: vector.hpp:150
Array< int > & GetSlave()
Definition: nurbs.hpp:320
void Generate3DBdrElementDofTable()
Definition: nurbs.cpp:2680
int GetNE() const
Returns number of elements.
Definition: mesh.hpp:676
int NumOfControlPoints
Definition: nurbs.hpp:37
Array< int > e_spaceOffsets
Definition: nurbs.hpp:200
void GenerateBdrElementDofTable()
Definition: nurbs.cpp:2618
const KnotVector * GetKnotVector(int i) const
Definition: nurbs.hpp:355
void GetElementLocalToGlobal(Array< int > &lelem_elem)
Definition: nurbs.cpp:2748
void init(int dim_)
Definition: nurbs.cpp:277
Array< int > edge_to_knot
Definition: nurbs.hpp:181
void GenerateActiveVertices()
Definition: nurbs.cpp:1846
Array< bool > activeElem
Definition: nurbs.hpp:175
void GetVertexLocalToGlobal(Array< int > &lvert_vert)
Definition: nurbs.cpp:2738
void Set2DSolutionVector(Vector &Nodes, int vdim)
Definition: nurbs.cpp:3066
Vector knot
Definition: nurbs.hpp:36
int GetNTotalDof() const
Definition: nurbs.hpp:351
friend NURBSPatch * Revolve3D(NURBSPatch &patch, double n[], double ang, int times)
Definition: nurbs.cpp:1120
int GetNDof() const
Definition: nurbs.hpp:352
int Dimension() const
Definition: nurbs.hpp:332
Array< int > el_to_patch
Definition: nurbs.hpp:206
Array< int > & GetMaster()
Definition: nurbs.hpp:318
void SwapDirections(int dir1, int dir2)
Definition: nurbs.cpp:960
KnotVector()
Create KnotVector.
Definition: nurbs.hpp:41
static const int MaxOrder
Definition: nurbs.hpp:34
Array< int > v_spaceOffsets
Definition: nurbs.hpp:199
void UniformRefinement(Vector &newknots) const
Definition: nurbs.cpp:86
const Array< int > & GetSlave() const
Definition: nurbs.hpp:319
void Generate2DBdrElementDofTable()
Definition: nurbs.cpp:2637
void ConnectBoundaries3D(int bnd0, int bnd1)
Definition: nurbs.cpp:1745
ParNURBSExtension(const ParNURBSExtension &orig)
Definition: nurbs.cpp:3127
Array< int > d_to_d
Definition: nurbs.hpp:188
void Set3DSolutionVector(Vector &Nodes, int vdim)
Definition: nurbs.cpp:3094
void GetPatchKnotVectors(int p, Array< KnotVector * > &kv)
Definition: nurbs.cpp:2059
void GenerateActiveBdrElems()
Definition: nurbs.cpp:1913
double * data
Definition: nurbs.hpp:89
void PrintSolution(const GridFunction &sol, std::ostream &out) const
Definition: nurbs.cpp:2896
void Get3DPatchNets(const Vector &Nodes, int vdim)
Definition: nurbs.cpp:3024
void mfem_error(const char *msg)
Function called when an error is encountered. Used by the macros MFEM_ABORT, MFEM_ASSERT, MFEM_VERIFY.
Definition: error.cpp:146
void Get2DPatchNets(const Vector &Nodes, int vdim)
Definition: nurbs.cpp:2997
Vector & GetWeights()
Definition: nurbs.hpp:374
bool HavePatches() const
Definition: nurbs.hpp:361
virtual ~NURBSExtension()
Destroy a NURBSExtension.
Definition: nurbs.cpp:1537
void GenerateOffsets()
Definition: nurbs.cpp:2163
Array< int > activeDof
Definition: nurbs.hpp:177
friend NURBSPatch * Interpolate(NURBSPatch &p1, NURBSPatch &p2)
Definition: nurbs.cpp:1073
void GetBdrElementTopo(Array< Element * > &boundary) const
Definition: nurbs.cpp:2385
static void Get3DRotationMatrix(double n[], double angle, double r, DenseMatrix &T)
Definition: nurbs.cpp:985
void SetKnotsFromPatches()
Definition: nurbs.cpp:2819
void GetElementTopo(Array< Element * > &elements) const
Definition: nurbs.cpp:2290
int GetNKV() const
Definition: nurbs.hpp:128
const Array< int > & GetOrders() const
Read-only access to the orders of all knot vectors.
Definition: nurbs.hpp:337
void UniformRefinement()
Definition: nurbs.cpp:541
Array< int > f_spaceOffsets
Definition: nurbs.hpp:201
int GetNBP() const
Definition: nurbs.hpp:334
Array< bool > activeBdrElem
Definition: nurbs.hpp:176
int GetNC() const
Definition: nurbs.hpp:127
void LoadSolution(std::istream &input, GridFunction &sol) const
Definition: nurbs.cpp:2859
int Dimension() const
Definition: mesh.hpp:713
KnotVector(const KnotVector &kv)
Definition: nurbs.hpp:44
int GetOrder() const
If all orders are identical, return that number. Otherwise, return NURBSFECollection::VariableOrder.
Definition: nurbs.hpp:340
KnotVector * GetKV(int i)
Definition: nurbs.hpp:129
Array< int > slave
Definition: nurbs.hpp:190
int GetGNV() const
Definition: nurbs.hpp:344
void Print(std::ostream &out) const
Definition: nurbs.cpp:458
void LoadFE(int i, const FiniteElement *FE) const
Definition: nurbs.cpp:2758
int DofMap(int dof) const
Definition: nurbs.hpp:234
Array< int > v_meshOffsets
Definition: nurbs.hpp:193
int GetNP() const
Definition: nurbs.hpp:333
void SetBdrPatchVertexMap(int p, const KnotVector *kv[], int *okv)
Definition: nurbs.cpp:3587
void GetBdrPatchKnotVectors(int p, Array< KnotVector * > &kv)
Definition: nurbs.cpp:2100
Array< int > p_meshOffsets
Definition: nurbs.hpp:196
KnotVector & operator=(const KnotVector &kv)
Definition: nurbs.cpp:45
void SetPatchVertexMap(int p, const KnotVector *kv[])
Definition: nurbs.cpp:3527
Table * GetElementDofTable()
Definition: nurbs.hpp:363
Table * GetBdrElementDofTable()
Definition: nurbs.hpp:364
Array2D< int > el_to_IJK
Definition: nurbs.hpp:208
void SetOrderFromOrders()
Definition: nurbs.cpp:2139
void CalcShape(Vector &shape, int i, double xi) const
Definition: nurbs.cpp:132
void CountBdrElements()
Definition: nurbs.cpp:2270
int GetNE() const
Definition: nurbs.hpp:347
void GetPatchNets(const Vector &Nodes, int vdim)
Definition: nurbs.cpp:2985
int MakeUniformDegree(int degree=-1)
Definition: nurbs.cpp:1050
bool isElement(int i) const
Definition: nurbs.hpp:57
Array< int > p_spaceOffsets
Definition: nurbs.hpp:202
void SetSolutionVector(Vector &Nodes, int vdim)
Definition: nurbs.cpp:3054
void Get3DElementTopo(Array< Element * > &elements) const
Definition: nurbs.cpp:2340
int findKnotSpan(double u) const
Definition: nurbs.cpp:214
void SetBdrPatchDofMap(int p, const KnotVector *kv[], int *okv)
Definition: nurbs.cpp:3616
void ConnectBoundaries()
Definition: nurbs.cpp:1640
int GetNCP() const
Definition: nurbs.hpp:50
double & operator()(int i, int j)
Definition: nurbs.hpp:474
void GenerateElementDofTable()
Definition: nurbs.cpp:2468
Array< NURBSPatch * > patches
Definition: nurbs.hpp:211
const double & operator[](int i) const
Definition: nurbs.hpp:81
void PrintCharacteristics(std::ostream &out) const
Definition: nurbs.cpp:1596
~KnotVector()
Destroys KnotVector.
Definition: nurbs.hpp:78
void Difference(const KnotVector &kv, Vector &diff) const
Definition: nurbs.cpp:243
Array< int > mOrders
Definition: nurbs.hpp:166
Vector data type.
Definition: vector.hpp:48
int operator[](const int i) const
Definition: nurbs.hpp:464
void Print(std::ostream &out) const
Definition: nurbs.cpp:125
double getKnotLocation(double xi, int ni) const
Definition: nurbs.hpp:59
int GetNBE() const
Definition: nurbs.hpp:349
int GetGNBE() const
Definition: nurbs.hpp:348
void CalcDShape(Vector &grad, int i, double xi) const
Definition: nurbs.cpp:158
void DegreeElevate(int rel_degree, int degree=16)
Definition: nurbs.cpp:2934
OutStream out(std::cout)
Global stream used by the library for standard output. Initially it uses the same std::streambuf as s...
Definition: globals.hpp:64
void CheckBdrPatches()
Definition: nurbs.cpp:2031
double & operator[](int i)
Definition: nurbs.hpp:80
void FlipDirection(int dir)
Definition: nurbs.cpp:948
void Get2DElementTopo(Array< Element * > &elements) const
Definition: nurbs.cpp:2304
GroupTopology gtopo
Definition: nurbs.hpp:409
int GetOrder() const
Definition: nurbs.hpp:51
Array< int > bel_to_patch
Definition: nurbs.hpp:207
void SetPatchDofMap(int p, const KnotVector *kv[])
Definition: nurbs.cpp:3557
void ConvertToPatches(const Vector &Nodes)
Definition: nurbs.cpp:2800
void KnotInsert(Array< KnotVector * > &kv)
Definition: nurbs.cpp:2958
void Print(std::ostream &out) const
Definition: nurbs.cpp:1561
void Get2DBdrElementTopo(Array< Element * > &boundary) const
Definition: nurbs.cpp:2399
int SetLoopDirection(int dir)
Definition: nurbs.cpp:482