MFEM  v4.5.1
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-2022, Lawrence Livermore National Security, LLC. Produced
2 // at the Lawrence Livermore National Laboratory. All Rights reserved. See files
3 // LICENSE and NOTICE for details. LLNL-CODE-806117.
4 //
5 // This file is part of the MFEM library. For more information and source code
6 // availability visit https://mfem.org.
7 //
8 // MFEM is free software; you can redistribute it and/or modify it under the
9 // terms of the BSD-3 license. We welcome feedback and contributions, see file
10 // CONTRIBUTING.md for details.
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  void CalcDnShape(Vector &gradn, int n, int i, double xi) const;
67  void CalcD2Shape(Vector &grad2, int i, double xi) const
68  { CalcDnShape(grad2, 2, i, xi); }
69 
70  void Difference(const KnotVector &kv, Vector &diff) const;
71  void UniformRefinement(Vector &newknots) const;
72  /** Return a new KnotVector with elevated degree by repeating the endpoints
73  of the knot vector. */
74  KnotVector *DegreeElevate(int t) const;
75 
76  void Flip();
77 
78  void Print(std::ostream &out) const;
79 
80  void PrintFunctions(std::ostream &out, int samples=11) const;
81 
82  /// Destroys KnotVector
84 
85  double &operator[](int i) { return knot(i); }
86  const double &operator[](int i) const { return knot(i); }
87 };
88 
89 
91 {
92 protected:
93  int ni, nj, nk, Dim;
94  double *data; // the layout of data is: (Dim x ni x nj x nk)
95 
97 
98  int sd, nd;
99 
100  void swap(NURBSPatch *np);
101 
102  // Special B-NET access functions
103  int SetLoopDirection(int dir);
104  inline double &operator()(int i, int j);
105  inline const double &operator()(int i, int j) const;
106 
107  void init(int dim_);
108 
109  NURBSPatch(NURBSPatch *parent, int dir, int Order, int NCP);
110 
111 public:
112  NURBSPatch(const NURBSPatch &orig);
113  NURBSPatch(std::istream &input);
114  NURBSPatch(const KnotVector *kv0, const KnotVector *kv1, int dim_);
115  NURBSPatch(const KnotVector *kv0, const KnotVector *kv1,
116  const KnotVector *kv2, int dim_);
118 
119  NURBSPatch& operator=(const NURBSPatch&) = delete;
120 
121  ~NURBSPatch();
122 
123  void Print(std::ostream &out) const;
124 
125  void DegreeElevate(int dir, int t);
126  void KnotInsert (int dir, const KnotVector &knot);
127  void KnotInsert (int dir, const Vector &knot);
128 
129  void KnotInsert(Array<Vector *> &knot);
130  void KnotInsert(Array<KnotVector *> &knot);
131 
132  void DegreeElevate(int t);
133  void UniformRefinement();
134 
135  // Return the number of components stored in the NURBSPatch
136  int GetNC() const { return Dim; }
137  int GetNKV() const { return kv.Size(); }
138  KnotVector *GetKV(int i) { return kv[i]; }
139 
140  // Standard B-NET access functions
141  inline double &operator()(int i, int j, int l);
142  inline const double &operator()(int i, int j, int l) const;
143 
144  inline double &operator()(int i, int j, int k, int l);
145  inline const double &operator()(int i, int j, int k, int l) const;
146 
147  static void Get3DRotationMatrix(double n[], double angle, double r,
148  DenseMatrix &T);
149  void FlipDirection(int dir);
150  void SwapDirections(int dir1, int dir2);
151  void Rotate3D(double normal[], double angle);
152  int MakeUniformDegree(int degree = -1);
153  friend NURBSPatch *Interpolate(NURBSPatch &p1, NURBSPatch &p2);
154  friend NURBSPatch *Revolve3D(NURBSPatch &patch, double n[], double ang,
155  int times);
156 };
157 
158 
159 #ifdef MFEM_USE_MPI
160 class ParNURBSExtension;
161 #endif
162 
163 class NURBSPatchMap;
164 
165 
167 {
168 #ifdef MFEM_USE_MPI
169  friend class ParNURBSExtension;
170 #endif
171  friend class NURBSPatchMap;
172 
173 protected:
174  int mOrder; // see GetOrder() for description
177  // global entity counts
179  // local entity counts
182 
183  Array<int> activeVert; // activeVert[glob_vert] = loc_vert or -1
186  Array<int> activeDof; // activeDof[glob_dof] = loc_dof + 1 or 0
187 
189  int own_topo;
193 
194  // periodic BC info:
195  // - dof 2 dof map
196  // - master and slave boundary indices
200 
201  // global offsets, meshOffsets == meshVertexOffsets
206 
207  // global offsets, spaceOffsets == dofOffsets
212 
214 
217  Array2D<int> el_to_IJK; // IJK are "knot-span" indices!
218  Array2D<int> bel_to_IJK; // they are NOT element indices!
219 
221 
222  inline int KnotInd(int edge) const;
223  inline KnotVector *KnotVec(int edge);
224  inline const KnotVector *KnotVec(int edge) const;
225  inline const KnotVector *KnotVec(int edge, int oedge, int *okv) const;
226 
227  void CheckPatches();
228  void CheckBdrPatches();
229 
230  void GetPatchKnotVectors (int p, Array<KnotVector *> &kv);
231  void GetPatchKnotVectors (int p, Array<const KnotVector *> &kv) const;
234 
235  void SetOrderFromOrders();
237 
238  // periodic BC helper functions
239  void InitDofMap();
240  void ConnectBoundaries();
241  void ConnectBoundaries2D(int bnd0, int bnd1);
242  void ConnectBoundaries3D(int bnd0, int bnd1);
243  int DofMap(int dof) const
244  {
245  return (d_to_d.Size() > 0 )? d_to_d[dof] : dof;
246  };
247 
248  // also count the global NumOfVertices and the global NumOfDofs
249  void GenerateOffsets();
250  // count the global NumOfElements
251  void CountElements();
252  // count the global NumOfBdrElements
253  void CountBdrElements();
254 
255  // generate the mesh elements
256  void Get2DElementTopo(Array<Element *> &elements) const;
257  void Get3DElementTopo(Array<Element *> &elements) const;
258 
259  // generate the boundary mesh elements
260  void Get2DBdrElementTopo(Array<Element *> &boundary) const;
261  void Get3DBdrElementTopo(Array<Element *> &boundary) const;
262 
263 
264  // FE space generation functions
265 
266  // based on activeElem, count NumOfActiveDofs, generate el_dof,
267  // el_to_patch, el_to_IJK, activeDof map (global-to-local)
269 
270  // generate elem_to_global-dof table for the active elements
271  // define el_to_patch, el_to_IJK, activeDof (as bool)
274 
275  // call after GenerateElementDofTable
277 
278  // generate the bdr-elem_to_global-dof table for the active bdr. elements
279  // define bel_to_patch, bel_to_IJK
282 
283  // FE --> Patch translation functions
284  void GetPatchNets (const Vector &Nodes, int vdim);
285  void Get2DPatchNets(const Vector &Nodes, int vdim);
286  void Get3DPatchNets(const Vector &Nodes, int vdim);
287 
288  // Patch --> FE translation functions
289  // Side effects: delete the patches, update the weights from the patches
290  void SetSolutionVector (Vector &Nodes, int vdim);
291  void Set2DSolutionVector(Vector &Nodes, int vdim);
292  void Set3DSolutionVector(Vector &Nodes, int vdim);
293 
294  // determine activeVert, NumOfActiveVertices from the activeElem array
295  void GenerateActiveVertices();
296 
297  // determine activeBdrElem, NumOfActiveBdrElems
298  void GenerateActiveBdrElems();
299 
300  void MergeWeights(Mesh *mesh_array[], int num_pieces);
301 
302  // to be used by ParNURBSExtension constructor(s)
304 
305 public:
306  /// Copy constructor: deep copy
307  NURBSExtension(const NURBSExtension &orig);
308  /// Read-in a NURBSExtension
309  NURBSExtension(std::istream &input);
310  /** @brief Create a NURBSExtension with elevated order by repeating the
311  endpoints of the knot vectors and using uniform weights of 1. */
312  /** If a knot vector in @a parent already has order greater than or equal to
313  @a newOrder, it will be used unmodified. */
314  NURBSExtension(NURBSExtension *parent, int newOrder);
315  /** @brief Create a NURBSExtension with elevated knot vector orders (by
316  repeating the endpoints of the knot vectors and using uniform weights of
317  1) as given by the array @a newOrders. */
318  /** If a knot vector in @a parent already has order greater than or equal to
319  the corresponding entry in @a newOrder, it will be used unmodified. */
320  NURBSExtension(NURBSExtension *parent, const Array<int> &newOrders);
321  /// Construct a NURBSExtension by merging a partitioned NURBS mesh
322  NURBSExtension(Mesh *mesh_array[], int num_pieces);
323 
324  /// Copy assignment not supported
325  NURBSExtension& operator=(const NURBSExtension&) = delete;
326 
327  // Generate connections between boundaries, such as periodic BCs
329  const Array<int> &GetMaster() const { return master; };
330  Array<int> &GetMaster() { return master; };
331  const Array<int> &GetSlave() const { return slave; };
332  Array<int> &GetSlave() { return slave; };
333  void MergeGridFunctions(GridFunction *gf_array[], int num_pieces,
334  GridFunction &merged);
335 
336  /// Destroy a NURBSExtension
337  virtual ~NURBSExtension();
338 
339  // Print functions
340  void Print(std::ostream &out) const;
341  void PrintCharacteristics(std::ostream &out) const;
342  void PrintFunctions(const char *filename, int samples=11) const;
343 
344  // Meta data functions
345  int Dimension() const { return patchTopo->Dimension(); }
346  int GetNP() const { return patchTopo->GetNE(); }
347  int GetNBP() const { return patchTopo->GetNBE(); }
348 
349  /// Read-only access to the orders of all knot vectors.
350  const Array<int> &GetOrders() const { return mOrders; }
351  /** @brief If all orders are identical, return that number. Otherwise, return
352  NURBSFECollection::VariableOrder. */
353  int GetOrder() const { return mOrder; }
354 
355  int GetNKV() const { return NumOfKnotVectors; }
356 
357  int GetGNV() const { return NumOfVertices; }
358  int GetNV() const { return NumOfActiveVertices; }
359  int GetGNE() const { return NumOfElements; }
360  int GetNE() const { return NumOfActiveElems; }
361  int GetGNBE() const { return NumOfBdrElements; }
362  int GetNBE() const { return NumOfActiveBdrElems; }
363 
364  int GetNTotalDof() const { return NumOfDofs; }
365  int GetNDof() const { return NumOfActiveDofs; }
366 
367  // Knotvector read-only access function
368  const KnotVector *GetKnotVector(int i) const { return knotVectors[i]; }
369 
370  // Mesh generation functions
371  void GetElementTopo (Array<Element *> &elements) const;
372  void GetBdrElementTopo(Array<Element *> &boundary) const;
373 
374  bool HavePatches() const { return (patches.Size() != 0); }
375 
378 
379  void GetVertexLocalToGlobal(Array<int> &lvert_vert);
380  void GetElementLocalToGlobal(Array<int> &lelem_elem);
381 
382  // Load functions
383  void LoadFE(int i, const FiniteElement *FE) const;
384  void LoadBE(int i, const FiniteElement *BE) const;
385 
386  const Vector &GetWeights() const { return weights; }
387  Vector &GetWeights() { return weights; }
388 
389  // Translation functions: from FE coordinates to IJK patch
390  // format and vice versa
391  void ConvertToPatches(const Vector &Nodes);
392  void SetKnotsFromPatches();
393  void SetCoordsFromPatches(Vector &Nodes);
394 
395  // Read a GridFunction written patch-by-patch, e.g. with PrintSolution().
396  void LoadSolution(std::istream &input, GridFunction &sol) const;
397  // Write a GridFunction patch-by-patch.
398  void PrintSolution(const GridFunction &sol, std::ostream &out) const;
399 
400  // Refinement methods
401  // new_degree = max(old_degree, min(old_degree + rel_degree, degree))
402  void DegreeElevate(int rel_degree, int degree = 16);
403  void UniformRefinement();
405  void KnotInsert(Array<Vector *> &kv);
406 };
407 
408 
409 #ifdef MFEM_USE_MPI
411 {
412 private:
413  int *partitioning;
414 
415  Table *GetGlobalElementDofTable();
416  Table *Get2DGlobalElementDofTable();
417  Table *Get3DGlobalElementDofTable();
418 
419  void SetActive(const int *partitioning, const Array<bool> &active_bel);
420  void BuildGroups(const int *partitioning, const Table &elem_dof);
421 
422 public:
424 
426 
428 
429  ParNURBSExtension(MPI_Comm comm, NURBSExtension *parent, int *partitioning,
430  const Array<bool> &active_bel);
431 
432  // Create a parallel version of 'parent' with partitioning as in
433  // 'par_parent'; the 'parent' object is destroyed.
434  // The 'parent' can be either a local NURBSExtension or a global one.
436  const ParNURBSExtension *par_parent);
437 
438  virtual ~ParNURBSExtension() { delete [] partitioning; }
439 };
440 #endif
441 
442 
444 {
445 private:
446  const NURBSExtension *Ext;
447 
448  int I, J, K, pOffset, opatch;
449  Array<int> verts, edges, faces, oedge, oface;
450 
451  inline static int F(const int n, const int N)
452  { return (n < 0) ? 0 : ((n >= N) ? 2 : 1); }
453 
454  inline static int Or1D(const int n, const int N, const int Or)
455  { return (Or > 0) ? n : (N - 1 - n); }
456 
457  inline static int Or2D(const int n1, const int n2,
458  const int N1, const int N2, const int Or);
459 
460  // also set verts, edges, faces, orientations etc
461  void GetPatchKnotVectors (int p, const KnotVector *kv[]);
462  void GetBdrPatchKnotVectors(int p, const KnotVector *kv[], int *okv);
463 
464 public:
465  NURBSPatchMap(const NURBSExtension *ext) { Ext = ext; }
466 
467  int nx() { return I + 1; }
468  int ny() { return J + 1; }
469  int nz() { return K + 1; }
470 
471  void SetPatchVertexMap(int p, const KnotVector *kv[]);
472  void SetPatchDofMap (int p, const KnotVector *kv[]);
473 
474  void SetBdrPatchVertexMap(int p, const KnotVector *kv[], int *okv);
475  void SetBdrPatchDofMap (int p, const KnotVector *kv[], int *okv);
476 
477  inline int operator()(const int i) const;
478  inline int operator[](const int i) const { return (*this)(i); }
479 
480  inline int operator()(const int i, const int j) const;
481 
482  inline int operator()(const int i, const int j, const int k) const;
483 };
484 
485 
486 // Inline function implementations
487 
488 inline double &NURBSPatch::operator()(int i, int j)
489 {
490  return data[j%sd + sd*(i + (j/sd)*nd)];
491 }
492 
493 inline const double &NURBSPatch::operator()(int i, int j) const
494 {
495  return data[j%sd + sd*(i + (j/sd)*nd)];
496 }
497 
498 inline double &NURBSPatch::operator()(int i, int j, int l)
499 {
500 #ifdef MFEM_DEBUG
501  if (data == 0 || i < 0 || i >= ni || j < 0 || j >= nj || nk > 0 ||
502  l < 0 || l >= Dim)
503  {
504  mfem_error("NURBSPatch::operator() 2D");
505  }
506 #endif
507 
508  return data[(i+j*ni)*Dim+l];
509 }
510 
511 inline const double &NURBSPatch::operator()(int i, int j, int l) const
512 {
513 #ifdef MFEM_DEBUG
514  if (data == 0 || i < 0 || i >= ni || j < 0 || j >= nj || nk > 0 ||
515  l < 0 || l >= Dim)
516  {
517  mfem_error("NURBSPatch::operator() const 2D");
518  }
519 #endif
520 
521  return data[(i+j*ni)*Dim+l];
522 }
523 
524 inline double &NURBSPatch::operator()(int i, int j, int k, int l)
525 {
526 #ifdef MFEM_DEBUG
527  if (data == 0 || i < 0 || i >= ni || j < 0 || j >= nj || k < 0 ||
528  k >= nk || l < 0 || l >= Dim)
529  {
530  mfem_error("NURBSPatch::operator() 3D");
531  }
532 #endif
533 
534  return data[(i+(j+k*nj)*ni)*Dim+l];
535 }
536 
537 inline const double &NURBSPatch::operator()(int i, int j, int k, int l) const
538 {
539 #ifdef MFEM_DEBUG
540  if (data == 0 || i < 0 || i >= ni || j < 0 || j >= nj || k < 0 ||
541  k >= nk || l < 0 || l >= Dim)
542  {
543  mfem_error("NURBSPatch::operator() const 3D");
544  }
545 #endif
546 
547  return data[(i+(j+k*nj)*ni)*Dim+l];
548 }
549 
550 
551 inline int NURBSExtension::KnotInd(int edge) const
552 {
553  int kv = edge_to_knot[edge];
554  return (kv >= 0) ? kv : (-1-kv);
555 }
556 
558 {
559  return knotVectors[KnotInd(edge)];
560 }
561 
562 inline const KnotVector *NURBSExtension::KnotVec(int edge) const
563 {
564  return knotVectors[KnotInd(edge)];
565 }
566 
567 inline const KnotVector *NURBSExtension::KnotVec(int edge, int oedge, int *okv)
568 const
569 {
570  int kv = edge_to_knot[edge];
571  if (kv >= 0)
572  {
573  *okv = oedge;
574  return knotVectors[kv];
575  }
576  else
577  {
578  *okv = -oedge;
579  return knotVectors[-1-kv];
580  }
581 }
582 
583 
584 // static method
585 inline int NURBSPatchMap::Or2D(const int n1, const int n2,
586  const int N1, const int N2, const int Or)
587 {
588  // Needs testing
589  switch (Or)
590  {
591  case 0: return n1 + n2*N1;
592  case 1: return n2 + n1*N2;
593  case 2: return n2 + (N1 - 1 - n1)*N2;
594  case 3: return (N1 - 1 - n1) + n2*N1;
595  case 4: return (N1 - 1 - n1) + (N2 - 1 - n2)*N1;
596  case 5: return (N2 - 1 - n2) + (N1 - 1 - n1)*N2;
597  case 6: return (N2 - 1 - n2) + n1*N2;
598  case 7: return n1 + (N2 - 1 - n2)*N1;
599  }
600 #ifdef MFEM_DEBUG
601  mfem_error("NURBSPatchMap::Or2D");
602 #endif
603  return -1;
604 }
605 
606 inline int NURBSPatchMap::operator()(const int i) const
607 {
608  const int i1 = i - 1;
609  switch (F(i1, I))
610  {
611  case 0: return verts[0];
612  case 1: return pOffset + Or1D(i1, I, opatch);
613  case 2: return verts[1];
614  }
615 #ifdef MFEM_DEBUG
616  mfem_error("NURBSPatchMap::operator() const 1D");
617 #endif
618  return -1;
619 }
620 
621 inline int NURBSPatchMap::operator()(const int i, const int j) const
622 {
623  const int i1 = i - 1, j1 = j - 1;
624  switch (3*F(j1, J) + F(i1, I))
625  {
626  case 0: return verts[0];
627  case 1: return edges[0] + Or1D(i1, I, oedge[0]);
628  case 2: return verts[1];
629  case 3: return edges[3] + Or1D(j1, J, -oedge[3]);
630  case 4: return pOffset + Or2D(i1, j1, I, J, opatch);
631  case 5: return edges[1] + Or1D(j1, J, oedge[1]);
632  case 6: return verts[3];
633  case 7: return edges[2] + Or1D(i1, I, -oedge[2]);
634  case 8: return verts[2];
635  }
636 #ifdef MFEM_DEBUG
637  mfem_error("NURBSPatchMap::operator() const 2D");
638 #endif
639  return -1;
640 }
641 
642 inline int NURBSPatchMap::operator()(const int i, const int j, const int k)
643 const
644 {
645  // Needs testing
646  const int i1 = i - 1, j1 = j - 1, k1 = k - 1;
647  switch (3*(3*F(k1, K) + F(j1, J)) + F(i1, I))
648  {
649  case 0: return verts[0];
650  case 1: return edges[0] + Or1D(i1, I, oedge[0]);
651  case 2: return verts[1];
652  case 3: return edges[3] + Or1D(j1, J, oedge[3]);
653  case 4: return faces[0] + Or2D(i1, J - 1 - j1, I, J, oface[0]);
654  case 5: return edges[1] + Or1D(j1, J, oedge[1]);
655  case 6: return verts[3];
656  case 7: return edges[2] + Or1D(i1, I, oedge[2]);
657  case 8: return verts[2];
658  case 9: return edges[8] + Or1D(k1, K, oedge[8]);
659  case 10: return faces[1] + Or2D(i1, k1, I, K, oface[1]);
660  case 11: return edges[9] + Or1D(k1, K, oedge[9]);
661  case 12: return faces[4] + Or2D(J - 1 - j1, k1, J, K, oface[4]);
662  case 13: return pOffset + I*(J*k1 + j1) + i1;
663  case 14: return faces[2] + Or2D(j1, k1, J, K, oface[2]);
664  case 15: return edges[11] + Or1D(k1, K, oedge[11]);
665  case 16: return faces[3] + Or2D(I - 1 - i1, k1, I, K, oface[3]);
666  case 17: return edges[10] + Or1D(k1, K, oedge[10]);
667  case 18: return verts[4];
668  case 19: return edges[4] + Or1D(i1, I, oedge[4]);
669  case 20: return verts[5];
670  case 21: return edges[7] + Or1D(j1, J, oedge[7]);
671  case 22: return faces[5] + Or2D(i1, j1, I, J, oface[5]);
672  case 23: return edges[5] + Or1D(j1, J, oedge[5]);
673  case 24: return verts[7];
674  case 25: return edges[6] + Or1D(i1, I, oedge[6]);
675  case 26: return verts[6];
676  }
677 #ifdef MFEM_DEBUG
678  mfem_error("NURBSPatchMap::operator() const 3D");
679 #endif
680  return -1;
681 }
682 
683 }
684 
685 #endif
Array< KnotVector * > knotVectors
Definition: nurbs.hpp:191
Abstract class for all finite elements.
Definition: fe_base.hpp:235
Array< KnotVector * > kv
Definition: nurbs.hpp:96
Array< int > f_meshOffsets
Definition: nurbs.hpp:204
KnotVector * DegreeElevate(int t) const
Definition: nurbs.cpp:58
NURBSPatch(NURBSPatch *parent, int dir, int Order, int NCP)
Definition: nurbs.cpp:534
Array2D< int > bel_to_IJK
Definition: nurbs.hpp:218
int Size() const
Return the logical size of the array.
Definition: array.hpp:138
void SetCoordsFromPatches(Vector &Nodes)
Definition: nurbs.cpp:2966
int Size() const
Definition: nurbs.hpp:52
int GetGNE() const
Definition: nurbs.hpp:359
Array< int > activeVert
Definition: nurbs.hpp:183
void ConnectBoundaries2D(int bnd0, int bnd1)
Definition: nurbs.cpp:1824
KnotVector * KnotVec(int edge)
Definition: nurbs.hpp:557
const Vector & GetWeights() const
Definition: nurbs.hpp:386
void swap(NURBSPatch *np)
Definition: nurbs.cpp:549
int GetNV() const
Definition: nurbs.hpp:358
const Array< int > & GetMaster() const
Definition: nurbs.hpp:329
Class for grid function - Vector with associated FE space.
Definition: gridfunc.hpp:30
void DegreeElevate(int dir, int t)
Definition: nurbs.cpp:837
int NumOfElements
Definition: nurbs.hpp:37
void LoadBE(int i, const FiniteElement *BE) const
Definition: nurbs.cpp:2934
void Generate2DElementDofTable()
Definition: nurbs.cpp:2653
int KnotInd(int edge) const
Definition: nurbs.hpp:551
void MergeGridFunctions(GridFunction *gf_array[], int num_pieces, GridFunction &merged)
Definition: nurbs.cpp:2115
Array< int > ldof_group
Definition: nurbs.hpp:425
Array< int > master
Definition: nurbs.hpp:198
virtual ~ParNURBSExtension()
Definition: nurbs.hpp:438
void Generate3DElementDofTable()
Definition: nurbs.cpp:2707
NURBSExtension & operator=(const NURBSExtension &)=delete
Copy assignment not supported.
void KnotInsert(int dir, const KnotVector &knot)
Definition: nurbs.cpp:689
void MergeWeights(Mesh *mesh_array[], int num_pieces)
Definition: nurbs.cpp:2090
int GetNBE() const
Returns number of boundary elements.
Definition: mesh.hpp:926
void Get3DBdrElementTopo(Array< Element * > &boundary) const
Definition: nurbs.cpp:2585
NURBSPatchMap(const NURBSExtension *ext)
Definition: nurbs.hpp:465
void UniformRefinement()
Definition: nurbs.cpp:3105
int GetNE() const
Definition: nurbs.hpp:48
void Rotate3D(double normal[], double angle)
Definition: nurbs.cpp:1164
Array< int > e_meshOffsets
Definition: nurbs.hpp:203
int GetNKV() const
Definition: nurbs.hpp:355
int GetNKS() const
Definition: nurbs.hpp:49
int operator()(const int i) const
Definition: nurbs.hpp:606
void GetElements()
Count the number of elements.
Definition: nurbs.cpp:101
void SetOrdersFromKnotVectors()
Definition: nurbs.cpp:2308
Data type dense matrix using column-major storage.
Definition: densemat.hpp:23
int Size() const
Returns the size of the vector.
Definition: vector.hpp:200
Array< int > & GetSlave()
Definition: nurbs.hpp:332
void Generate3DBdrElementDofTable()
Definition: nurbs.cpp:2835
int GetNE() const
Returns number of elements.
Definition: mesh.hpp:923
int NumOfControlPoints
Definition: nurbs.hpp:37
Array< int > e_spaceOffsets
Definition: nurbs.hpp:209
void GenerateBdrElementDofTable()
Definition: nurbs.cpp:2773
const KnotVector * GetKnotVector(int i) const
Definition: nurbs.hpp:368
void GetElementLocalToGlobal(Array< int > &lelem_elem)
Definition: nurbs.cpp:2903
void init(int dim_)
Definition: nurbs.cpp:407
Array< int > edge_to_knot
Definition: nurbs.hpp:190
void GenerateActiveVertices()
Definition: nurbs.cpp:2001
Array< bool > activeElem
Definition: nurbs.hpp:184
void GetVertexLocalToGlobal(Array< int > &lvert_vert)
Definition: nurbs.cpp:2893
void Set2DSolutionVector(Vector &Nodes, int vdim)
Definition: nurbs.cpp:3249
Vector knot
Definition: nurbs.hpp:36
int GetNTotalDof() const
Definition: nurbs.hpp:364
friend NURBSPatch * Revolve3D(NURBSPatch &patch, double n[], double ang, int times)
Definition: nurbs.cpp:1260
int GetNDof() const
Definition: nurbs.hpp:365
int Dimension() const
Definition: nurbs.hpp:345
Array< int > el_to_patch
Definition: nurbs.hpp:215
Array< int > & GetMaster()
Definition: nurbs.hpp:330
void SwapDirections(int dir1, int dir2)
Definition: nurbs.cpp:1100
KnotVector()
Create KnotVector.
Definition: nurbs.hpp:41
static const int MaxOrder
Definition: nurbs.hpp:34
Array< int > v_spaceOffsets
Definition: nurbs.hpp:208
void UniformRefinement(Vector &newknots) const
Definition: nurbs.cpp:87
const Array< int > & GetSlave() const
Definition: nurbs.hpp:331
void Generate2DBdrElementDofTable()
Definition: nurbs.cpp:2792
void ConnectBoundaries3D(int bnd0, int bnd1)
Definition: nurbs.cpp:1900
ParNURBSExtension(const ParNURBSExtension &orig)
Definition: nurbs.cpp:3310
Array< int > d_to_d
Definition: nurbs.hpp:197
void Set3DSolutionVector(Vector &Nodes, int vdim)
Definition: nurbs.cpp:3277
void GetPatchKnotVectors(int p, Array< KnotVector * > &kv)
Definition: nurbs.cpp:2214
void GenerateActiveBdrElems()
Definition: nurbs.cpp:2068
double * data
Definition: nurbs.hpp:94
void PrintSolution(const GridFunction &sol, std::ostream &out) const
Definition: nurbs.cpp:3051
void Get3DPatchNets(const Vector &Nodes, int vdim)
Definition: nurbs.cpp:3207
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:154
void Get2DPatchNets(const Vector &Nodes, int vdim)
Definition: nurbs.cpp:3180
Vector & GetWeights()
Definition: nurbs.hpp:387
bool HavePatches() const
Definition: nurbs.hpp:374
virtual ~NURBSExtension()
Destroy a NURBSExtension.
Definition: nurbs.cpp:1677
void GenerateOffsets()
Definition: nurbs.cpp:2318
Array< int > activeDof
Definition: nurbs.hpp:186
friend NURBSPatch * Interpolate(NURBSPatch &p1, NURBSPatch &p2)
Definition: nurbs.cpp:1213
void GetBdrElementTopo(Array< Element * > &boundary) const
Definition: nurbs.cpp:2540
static void Get3DRotationMatrix(double n[], double angle, double r, DenseMatrix &T)
Definition: nurbs.cpp:1125
void SetKnotsFromPatches()
Definition: nurbs.cpp:2974
void GetElementTopo(Array< Element * > &elements) const
Definition: nurbs.cpp:2445
int GetNKV() const
Definition: nurbs.hpp:137
const Array< int > & GetOrders() const
Read-only access to the orders of all knot vectors.
Definition: nurbs.hpp:350
void UniformRefinement()
Definition: nurbs.cpp:671
Array< int > f_spaceOffsets
Definition: nurbs.hpp:210
int GetNBP() const
Definition: nurbs.hpp:347
Array< bool > activeBdrElem
Definition: nurbs.hpp:185
int GetNC() const
Definition: nurbs.hpp:136
void LoadSolution(std::istream &input, GridFunction &sol) const
Definition: nurbs.cpp:3014
int Dimension() const
Definition: mesh.hpp:1006
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:353
KnotVector * GetKV(int i)
Definition: nurbs.hpp:138
Array< int > slave
Definition: nurbs.hpp:199
int GetGNV() const
Definition: nurbs.hpp:357
void Print(std::ostream &out) const
Definition: nurbs.cpp:588
void LoadFE(int i, const FiniteElement *FE) const
Definition: nurbs.cpp:2913
int DofMap(int dof) const
Definition: nurbs.hpp:243
Array< int > v_meshOffsets
Definition: nurbs.hpp:202
int GetNP() const
Definition: nurbs.hpp:346
void SetBdrPatchVertexMap(int p, const KnotVector *kv[], int *okv)
Definition: nurbs.cpp:3770
void GetBdrPatchKnotVectors(int p, Array< KnotVector * > &kv)
Definition: nurbs.cpp:2255
Array< int > p_meshOffsets
Definition: nurbs.hpp:205
KnotVector & operator=(const KnotVector &kv)
Definition: nurbs.cpp:46
void SetPatchVertexMap(int p, const KnotVector *kv[])
Definition: nurbs.cpp:3710
Table * GetElementDofTable()
Definition: nurbs.hpp:376
Table * GetBdrElementDofTable()
Definition: nurbs.hpp:377
Array2D< int > el_to_IJK
Definition: nurbs.hpp:217
void SetOrderFromOrders()
Definition: nurbs.cpp:2294
void CalcD2Shape(Vector &grad2, int i, double xi) const
Definition: nurbs.hpp:67
void CalcShape(Vector &shape, int i, double xi) const
Definition: nurbs.cpp:160
void CalcDnShape(Vector &gradn, int n, int i, double xi) const
Definition: nurbs.cpp:243
void CountBdrElements()
Definition: nurbs.cpp:2425
void PrintFunctions(std::ostream &out, int samples=11) const
Definition: nurbs.cpp:133
int GetNE() const
Definition: nurbs.hpp:360
void GetPatchNets(const Vector &Nodes, int vdim)
Definition: nurbs.cpp:3168
int MakeUniformDegree(int degree=-1)
Definition: nurbs.cpp:1190
bool isElement(int i) const
Definition: nurbs.hpp:57
Array< int > p_spaceOffsets
Definition: nurbs.hpp:211
void SetSolutionVector(Vector &Nodes, int vdim)
Definition: nurbs.cpp:3237
void Get3DElementTopo(Array< Element * > &elements) const
Definition: nurbs.cpp:2495
int findKnotSpan(double u) const
Definition: nurbs.cpp:345
void SetBdrPatchDofMap(int p, const KnotVector *kv[], int *okv)
Definition: nurbs.cpp:3799
void ConnectBoundaries()
Definition: nurbs.cpp:1793
int GetNCP() const
Definition: nurbs.hpp:50
double & operator()(int i, int j)
Definition: nurbs.hpp:488
void GenerateElementDofTable()
Definition: nurbs.cpp:2623
Array< NURBSPatch * > patches
Definition: nurbs.hpp:220
const double & operator[](int i) const
Definition: nurbs.hpp:86
void PrintCharacteristics(std::ostream &out) const
Definition: nurbs.cpp:1736
~KnotVector()
Destroys KnotVector.
Definition: nurbs.hpp:83
void Difference(const KnotVector &kv, Vector &diff) const
Definition: nurbs.cpp:374
NURBSPatch & operator=(const NURBSPatch &)=delete
Array< int > mOrders
Definition: nurbs.hpp:175
RefCoord t[3]
Vector data type.
Definition: vector.hpp:60
int operator[](const int i) const
Definition: nurbs.hpp:478
void Print(std::ostream &out) const
Definition: nurbs.cpp:126
double getKnotLocation(double xi, int ni) const
Definition: nurbs.hpp:59
int GetNBE() const
Definition: nurbs.hpp:362
void PrintFunctions(const char *filename, int samples=11) const
Definition: nurbs.cpp:1766
double u(const Vector &xvec)
Definition: lor_mms.hpp:24
int GetGNBE() const
Definition: nurbs.hpp:361
void CalcDShape(Vector &grad, int i, double xi) const
Definition: nurbs.cpp:186
void DegreeElevate(int rel_degree, int degree=16)
Definition: nurbs.cpp:3089
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:66
void CheckBdrPatches()
Definition: nurbs.cpp:2186
double & operator[](int i)
Definition: nurbs.hpp:85
void FlipDirection(int dir)
Definition: nurbs.cpp:1088
void Get2DElementTopo(Array< Element * > &elements) const
Definition: nurbs.cpp:2459
GroupTopology gtopo
Definition: nurbs.hpp:423
int GetOrder() const
Definition: nurbs.hpp:51
Array< int > bel_to_patch
Definition: nurbs.hpp:216
void SetPatchDofMap(int p, const KnotVector *kv[])
Definition: nurbs.cpp:3740
void ConvertToPatches(const Vector &Nodes)
Definition: nurbs.cpp:2955
void KnotInsert(Array< KnotVector * > &kv)
Definition: nurbs.cpp:3113
void Print(std::ostream &out) const
Definition: nurbs.cpp:1701
void Get2DBdrElementTopo(Array< Element * > &boundary) const
Definition: nurbs.cpp:2554
int SetLoopDirection(int dir)
Definition: nurbs.cpp:612