MFEM  v4.3.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-2021, 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();
120 
121  void Print(std::ostream &out) const;
122 
123  void DegreeElevate(int dir, int t);
124  void KnotInsert (int dir, const KnotVector &knot);
125  void KnotInsert (int dir, const Vector &knot);
126 
127  void KnotInsert(Array<Vector *> &knot);
128  void KnotInsert(Array<KnotVector *> &knot);
129 
130  void DegreeElevate(int t);
131  void UniformRefinement();
132 
133  // Return the number of components stored in the NURBSPatch
134  int GetNC() const { return Dim; }
135  int GetNKV() const { return kv.Size(); }
136  KnotVector *GetKV(int i) { return kv[i]; }
137 
138  // Standard B-NET access functions
139  inline double &operator()(int i, int j, int l);
140  inline const double &operator()(int i, int j, int l) const;
141 
142  inline double &operator()(int i, int j, int k, int l);
143  inline const double &operator()(int i, int j, int k, int l) const;
144 
145  static void Get3DRotationMatrix(double n[], double angle, double r,
146  DenseMatrix &T);
147  void FlipDirection(int dir);
148  void SwapDirections(int dir1, int dir2);
149  void Rotate3D(double normal[], double angle);
150  int MakeUniformDegree(int degree = -1);
151  friend NURBSPatch *Interpolate(NURBSPatch &p1, NURBSPatch &p2);
152  friend NURBSPatch *Revolve3D(NURBSPatch &patch, double n[], double ang,
153  int times);
154 };
155 
156 
157 #ifdef MFEM_USE_MPI
158 class ParNURBSExtension;
159 #endif
160 
161 class NURBSPatchMap;
162 
163 
165 {
166 #ifdef MFEM_USE_MPI
167  friend class ParNURBSExtension;
168 #endif
169  friend class NURBSPatchMap;
170 
171 protected:
172  int mOrder; // see GetOrder() for description
175  // global entity counts
177  // local entity counts
180 
181  Array<int> activeVert; // activeVert[glob_vert] = loc_vert or -1
184  Array<int> activeDof; // activeDof[glob_dof] = loc_dof + 1 or 0
185 
187  int own_topo;
191 
192  // periodic BC info:
193  // - dof 2 dof map
194  // - master and slave boundary indices
198 
199  // global offsets, meshOffsets == meshVertexOffsets
204 
205  // global offsets, spaceOffsets == dofOffsets
210 
212 
215  Array2D<int> el_to_IJK; // IJK are "knot-span" indices!
216  Array2D<int> bel_to_IJK; // they are NOT element indices!
217 
219 
220  inline int KnotInd(int edge) const;
221  inline KnotVector *KnotVec(int edge);
222  inline const KnotVector *KnotVec(int edge) const;
223  inline const KnotVector *KnotVec(int edge, int oedge, int *okv) const;
224 
225  void CheckPatches();
226  void CheckBdrPatches();
227 
228  void GetPatchKnotVectors (int p, Array<KnotVector *> &kv);
229  void GetPatchKnotVectors (int p, Array<const KnotVector *> &kv) const;
232 
233  void SetOrderFromOrders();
235 
236  // periodic BC helper functions
237  void InitDofMap();
238  void ConnectBoundaries();
239  void ConnectBoundaries2D(int bnd0, int bnd1);
240  void ConnectBoundaries3D(int bnd0, int bnd1);
241  int DofMap(int dof) const
242  {
243  return (d_to_d.Size() > 0 )? d_to_d[dof] : dof;
244  };
245 
246  // also count the global NumOfVertices and the global NumOfDofs
247  void GenerateOffsets();
248  // count the global NumOfElements
249  void CountElements();
250  // count the global NumOfBdrElements
251  void CountBdrElements();
252 
253  // generate the mesh elements
254  void Get2DElementTopo(Array<Element *> &elements) const;
255  void Get3DElementTopo(Array<Element *> &elements) const;
256 
257  // generate the boundary mesh elements
258  void Get2DBdrElementTopo(Array<Element *> &boundary) const;
259  void Get3DBdrElementTopo(Array<Element *> &boundary) const;
260 
261 
262  // FE space generation functions
263 
264  // based on activeElem, count NumOfActiveDofs, generate el_dof,
265  // el_to_patch, el_to_IJK, activeDof map (global-to-local)
267 
268  // generate elem_to_global-dof table for the active elements
269  // define el_to_patch, el_to_IJK, activeDof (as bool)
272 
273  // call after GenerateElementDofTable
275 
276  // generate the bdr-elem_to_global-dof table for the active bdr. elements
277  // define bel_to_patch, bel_to_IJK
280 
281  // FE --> Patch translation functions
282  void GetPatchNets (const Vector &Nodes, int vdim);
283  void Get2DPatchNets(const Vector &Nodes, int vdim);
284  void Get3DPatchNets(const Vector &Nodes, int vdim);
285 
286  // Patch --> FE translation functions
287  // Side effects: delete the patches, update the weights from the patches
288  void SetSolutionVector (Vector &Nodes, int vdim);
289  void Set2DSolutionVector(Vector &Nodes, int vdim);
290  void Set3DSolutionVector(Vector &Nodes, int vdim);
291 
292  // determine activeVert, NumOfActiveVertices from the activeElem array
293  void GenerateActiveVertices();
294 
295  // determine activeBdrElem, NumOfActiveBdrElems
296  void GenerateActiveBdrElems();
297 
298  void MergeWeights(Mesh *mesh_array[], int num_pieces);
299 
300  // to be used by ParNURBSExtension constructor(s)
302 
303 public:
304  /// Copy constructor: deep copy
305  NURBSExtension(const NURBSExtension &orig);
306  /// Read-in a NURBSExtension
307  NURBSExtension(std::istream &input);
308  /** @brief Create a NURBSExtension with elevated order by repeating the
309  endpoints of the knot vectors and using uniform weights of 1. */
310  /** If a knot vector in @a parent already has order greater than or equal to
311  @a newOrder, it will be used unmodified. */
312  NURBSExtension(NURBSExtension *parent, int newOrder);
313  /** @brief Create a NURBSExtension with elevated knot vector orders (by
314  repeating the endpoints of the knot vectors and using uniform weights of
315  1) as given by the array @a newOrders. */
316  /** If a knot vector in @a parent already has order greater than or equal to
317  the corresponding entry in @a newOrder, it will be used unmodified. */
318  NURBSExtension(NURBSExtension *parent, const Array<int> &newOrders);
319  /// Construct a NURBSExtension by merging a partitioned NURBS mesh
320  NURBSExtension(Mesh *mesh_array[], int num_pieces);
321 
322  // Generate connections between boundaries, such as periodic BCs
324  const Array<int> &GetMaster() const { return master; };
325  Array<int> &GetMaster() { return master; };
326  const Array<int> &GetSlave() const { return slave; };
327  Array<int> &GetSlave() { return slave; };
328  void MergeGridFunctions(GridFunction *gf_array[], int num_pieces,
329  GridFunction &merged);
330 
331  /// Destroy a NURBSExtension
332  virtual ~NURBSExtension();
333 
334  // Print functions
335  void Print(std::ostream &out) const;
336  void PrintCharacteristics(std::ostream &out) const;
337  void PrintFunctions(const char *filename, int samples=11) const;
338 
339  // Meta data functions
340  int Dimension() const { return patchTopo->Dimension(); }
341  int GetNP() const { return patchTopo->GetNE(); }
342  int GetNBP() const { return patchTopo->GetNBE(); }
343 
344  /// Read-only access to the orders of all knot vectors.
345  const Array<int> &GetOrders() const { return mOrders; }
346  /** @brief If all orders are identical, return that number. Otherwise, return
347  NURBSFECollection::VariableOrder. */
348  int GetOrder() const { return mOrder; }
349 
350  int GetNKV() const { return NumOfKnotVectors; }
351 
352  int GetGNV() const { return NumOfVertices; }
353  int GetNV() const { return NumOfActiveVertices; }
354  int GetGNE() const { return NumOfElements; }
355  int GetNE() const { return NumOfActiveElems; }
356  int GetGNBE() const { return NumOfBdrElements; }
357  int GetNBE() const { return NumOfActiveBdrElems; }
358 
359  int GetNTotalDof() const { return NumOfDofs; }
360  int GetNDof() const { return NumOfActiveDofs; }
361 
362  // Knotvector read-only access function
363  const KnotVector *GetKnotVector(int i) const { return knotVectors[i]; }
364 
365  // Mesh generation functions
366  void GetElementTopo (Array<Element *> &elements) const;
367  void GetBdrElementTopo(Array<Element *> &boundary) const;
368 
369  bool HavePatches() const { return (patches.Size() != 0); }
370 
373 
374  void GetVertexLocalToGlobal(Array<int> &lvert_vert);
375  void GetElementLocalToGlobal(Array<int> &lelem_elem);
376 
377  // Load functions
378  void LoadFE(int i, const FiniteElement *FE) const;
379  void LoadBE(int i, const FiniteElement *BE) const;
380 
381  const Vector &GetWeights() const { return weights; }
382  Vector &GetWeights() { return weights; }
383 
384  // Translation functions: from FE coordinates to IJK patch
385  // format and vice versa
386  void ConvertToPatches(const Vector &Nodes);
387  void SetKnotsFromPatches();
388  void SetCoordsFromPatches(Vector &Nodes);
389 
390  // Read a GridFunction written patch-by-patch, e.g. with PrintSolution().
391  void LoadSolution(std::istream &input, GridFunction &sol) const;
392  // Write a GridFunction patch-by-patch.
393  void PrintSolution(const GridFunction &sol, std::ostream &out) const;
394 
395  // Refinement methods
396  // new_degree = max(old_degree, min(old_degree + rel_degree, degree))
397  void DegreeElevate(int rel_degree, int degree = 16);
398  void UniformRefinement();
400  void KnotInsert(Array<Vector *> &kv);
401 };
402 
403 
404 #ifdef MFEM_USE_MPI
406 {
407 private:
408  int *partitioning;
409 
410  Table *GetGlobalElementDofTable();
411  Table *Get2DGlobalElementDofTable();
412  Table *Get3DGlobalElementDofTable();
413 
414  void SetActive(const int *partitioning, const Array<bool> &active_bel);
415  void BuildGroups(const int *partitioning, const Table &elem_dof);
416 
417 public:
419 
421 
423 
424  ParNURBSExtension(MPI_Comm comm, NURBSExtension *parent, int *partitioning,
425  const Array<bool> &active_bel);
426 
427  // Create a parallel version of 'parent' with partitioning as in
428  // 'par_parent'; the 'parent' object is destroyed.
429  // The 'parent' can be either a local NURBSExtension or a global one.
431  const ParNURBSExtension *par_parent);
432 
433  virtual ~ParNURBSExtension() { delete [] partitioning; }
434 };
435 #endif
436 
437 
439 {
440 private:
441  const NURBSExtension *Ext;
442 
443  int I, J, K, pOffset, opatch;
444  Array<int> verts, edges, faces, oedge, oface;
445 
446  inline static int F(const int n, const int N)
447  { return (n < 0) ? 0 : ((n >= N) ? 2 : 1); }
448 
449  inline static int Or1D(const int n, const int N, const int Or)
450  { return (Or > 0) ? n : (N - 1 - n); }
451 
452  inline static int Or2D(const int n1, const int n2,
453  const int N1, const int N2, const int Or);
454 
455  // also set verts, edges, faces, orientations etc
456  void GetPatchKnotVectors (int p, const KnotVector *kv[]);
457  void GetBdrPatchKnotVectors(int p, const KnotVector *kv[], int *okv);
458 
459 public:
460  NURBSPatchMap(const NURBSExtension *ext) { Ext = ext; }
461 
462  int nx() { return I + 1; }
463  int ny() { return J + 1; }
464  int nz() { return K + 1; }
465 
466  void SetPatchVertexMap(int p, const KnotVector *kv[]);
467  void SetPatchDofMap (int p, const KnotVector *kv[]);
468 
469  void SetBdrPatchVertexMap(int p, const KnotVector *kv[], int *okv);
470  void SetBdrPatchDofMap (int p, const KnotVector *kv[], int *okv);
471 
472  inline int operator()(const int i) const;
473  inline int operator[](const int i) const { return (*this)(i); }
474 
475  inline int operator()(const int i, const int j) const;
476 
477  inline int operator()(const int i, const int j, const int k) const;
478 };
479 
480 
481 // Inline function implementations
482 
483 inline double &NURBSPatch::operator()(int i, int j)
484 {
485  return data[j%sd + sd*(i + (j/sd)*nd)];
486 }
487 
488 inline const double &NURBSPatch::operator()(int i, int j) const
489 {
490  return data[j%sd + sd*(i + (j/sd)*nd)];
491 }
492 
493 inline double &NURBSPatch::operator()(int i, int j, int l)
494 {
495 #ifdef MFEM_DEBUG
496  if (data == 0 || i < 0 || i >= ni || j < 0 || j >= nj || nk > 0 ||
497  l < 0 || l >= Dim)
498  {
499  mfem_error("NURBSPatch::operator() 2D");
500  }
501 #endif
502 
503  return data[(i+j*ni)*Dim+l];
504 }
505 
506 inline const double &NURBSPatch::operator()(int i, int j, int l) const
507 {
508 #ifdef MFEM_DEBUG
509  if (data == 0 || i < 0 || i >= ni || j < 0 || j >= nj || nk > 0 ||
510  l < 0 || l >= Dim)
511  {
512  mfem_error("NURBSPatch::operator() const 2D");
513  }
514 #endif
515 
516  return data[(i+j*ni)*Dim+l];
517 }
518 
519 inline double &NURBSPatch::operator()(int i, int j, int k, int l)
520 {
521 #ifdef MFEM_DEBUG
522  if (data == 0 || i < 0 || i >= ni || j < 0 || j >= nj || k < 0 ||
523  k >= nk || l < 0 || l >= Dim)
524  {
525  mfem_error("NURBSPatch::operator() 3D");
526  }
527 #endif
528 
529  return data[(i+(j+k*nj)*ni)*Dim+l];
530 }
531 
532 inline const double &NURBSPatch::operator()(int i, int j, int k, int l) const
533 {
534 #ifdef MFEM_DEBUG
535  if (data == 0 || i < 0 || i >= ni || j < 0 || j >= nj || k < 0 ||
536  k >= nk || l < 0 || l >= Dim)
537  {
538  mfem_error("NURBSPatch::operator() const 3D");
539  }
540 #endif
541 
542  return data[(i+(j+k*nj)*ni)*Dim+l];
543 }
544 
545 
546 inline int NURBSExtension::KnotInd(int edge) const
547 {
548  int kv = edge_to_knot[edge];
549  return (kv >= 0) ? kv : (-1-kv);
550 }
551 
553 {
554  return knotVectors[KnotInd(edge)];
555 }
556 
557 inline const KnotVector *NURBSExtension::KnotVec(int edge) const
558 {
559  return knotVectors[KnotInd(edge)];
560 }
561 
562 inline const KnotVector *NURBSExtension::KnotVec(int edge, int oedge, int *okv)
563 const
564 {
565  int kv = edge_to_knot[edge];
566  if (kv >= 0)
567  {
568  *okv = oedge;
569  return knotVectors[kv];
570  }
571  else
572  {
573  *okv = -oedge;
574  return knotVectors[-1-kv];
575  }
576 }
577 
578 
579 // static method
580 inline int NURBSPatchMap::Or2D(const int n1, const int n2,
581  const int N1, const int N2, const int Or)
582 {
583  // Needs testing
584  switch (Or)
585  {
586  case 0: return n1 + n2*N1;
587  case 1: return n2 + n1*N2;
588  case 2: return n2 + (N1 - 1 - n1)*N2;
589  case 3: return (N1 - 1 - n1) + n2*N1;
590  case 4: return (N1 - 1 - n1) + (N2 - 1 - n2)*N1;
591  case 5: return (N2 - 1 - n2) + (N1 - 1 - n1)*N2;
592  case 6: return (N2 - 1 - n2) + n1*N2;
593  case 7: return n1 + (N2 - 1 - n2)*N1;
594  }
595 #ifdef MFEM_DEBUG
596  mfem_error("NURBSPatchMap::Or2D");
597 #endif
598  return -1;
599 }
600 
601 inline int NURBSPatchMap::operator()(const int i) const
602 {
603  const int i1 = i - 1;
604  switch (F(i1, I))
605  {
606  case 0: return verts[0];
607  case 1: return pOffset + Or1D(i1, I, opatch);
608  case 2: return verts[1];
609  }
610 #ifdef MFEM_DEBUG
611  mfem_error("NURBSPatchMap::operator() const 1D");
612 #endif
613  return -1;
614 }
615 
616 inline int NURBSPatchMap::operator()(const int i, const int j) const
617 {
618  const int i1 = i - 1, j1 = j - 1;
619  switch (3*F(j1, J) + F(i1, I))
620  {
621  case 0: return verts[0];
622  case 1: return edges[0] + Or1D(i1, I, oedge[0]);
623  case 2: return verts[1];
624  case 3: return edges[3] + Or1D(j1, J, -oedge[3]);
625  case 4: return pOffset + Or2D(i1, j1, I, J, opatch);
626  case 5: return edges[1] + Or1D(j1, J, oedge[1]);
627  case 6: return verts[3];
628  case 7: return edges[2] + Or1D(i1, I, -oedge[2]);
629  case 8: return verts[2];
630  }
631 #ifdef MFEM_DEBUG
632  mfem_error("NURBSPatchMap::operator() const 2D");
633 #endif
634  return -1;
635 }
636 
637 inline int NURBSPatchMap::operator()(const int i, const int j, const int k)
638 const
639 {
640  // Needs testing
641  const int i1 = i - 1, j1 = j - 1, k1 = k - 1;
642  switch (3*(3*F(k1, K) + F(j1, J)) + F(i1, I))
643  {
644  case 0: return verts[0];
645  case 1: return edges[0] + Or1D(i1, I, oedge[0]);
646  case 2: return verts[1];
647  case 3: return edges[3] + Or1D(j1, J, oedge[3]);
648  case 4: return faces[0] + Or2D(i1, J - 1 - j1, I, J, oface[0]);
649  case 5: return edges[1] + Or1D(j1, J, oedge[1]);
650  case 6: return verts[3];
651  case 7: return edges[2] + Or1D(i1, I, oedge[2]);
652  case 8: return verts[2];
653  case 9: return edges[8] + Or1D(k1, K, oedge[8]);
654  case 10: return faces[1] + Or2D(i1, k1, I, K, oface[1]);
655  case 11: return edges[9] + Or1D(k1, K, oedge[9]);
656  case 12: return faces[4] + Or2D(J - 1 - j1, k1, J, K, oface[4]);
657  case 13: return pOffset + I*(J*k1 + j1) + i1;
658  case 14: return faces[2] + Or2D(j1, k1, J, K, oface[2]);
659  case 15: return edges[11] + Or1D(k1, K, oedge[11]);
660  case 16: return faces[3] + Or2D(I - 1 - i1, k1, I, K, oface[3]);
661  case 17: return edges[10] + Or1D(k1, K, oedge[10]);
662  case 18: return verts[4];
663  case 19: return edges[4] + Or1D(i1, I, oedge[4]);
664  case 20: return verts[5];
665  case 21: return edges[7] + Or1D(j1, J, oedge[7]);
666  case 22: return faces[5] + Or2D(i1, j1, I, J, oface[5]);
667  case 23: return edges[5] + Or1D(j1, J, oedge[5]);
668  case 24: return verts[7];
669  case 25: return edges[6] + Or1D(i1, I, oedge[6]);
670  case 26: return verts[6];
671  }
672 #ifdef MFEM_DEBUG
673  mfem_error("NURBSPatchMap::operator() const 3D");
674 #endif
675  return -1;
676 }
677 
678 }
679 
680 #endif
Array< KnotVector * > knotVectors
Definition: nurbs.hpp:189
Abstract class for all finite elements.
Definition: fe.hpp:243
Array< KnotVector * > kv
Definition: nurbs.hpp:96
Array< int > f_meshOffsets
Definition: nurbs.hpp:202
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:216
int Size() const
Return the logical size of the array.
Definition: array.hpp:134
void SetCoordsFromPatches(Vector &Nodes)
Definition: nurbs.cpp:2966
int Size() const
Definition: nurbs.hpp:52
int GetGNE() const
Definition: nurbs.hpp:354
Array< int > activeVert
Definition: nurbs.hpp:181
void ConnectBoundaries2D(int bnd0, int bnd1)
Definition: nurbs.cpp:1824
KnotVector * KnotVec(int edge)
Definition: nurbs.hpp:552
const Vector & GetWeights() const
Definition: nurbs.hpp:381
void swap(NURBSPatch *np)
Definition: nurbs.cpp:549
int GetNV() const
Definition: nurbs.hpp:353
const Array< int > & GetMaster() const
Definition: nurbs.hpp:324
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:546
void MergeGridFunctions(GridFunction *gf_array[], int num_pieces, GridFunction &merged)
Definition: nurbs.cpp:2115
Array< int > ldof_group
Definition: nurbs.hpp:420
Array< int > master
Definition: nurbs.hpp:196
virtual ~ParNURBSExtension()
Definition: nurbs.hpp:433
void Generate3DElementDofTable()
Definition: nurbs.cpp:2707
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:849
void Get3DBdrElementTopo(Array< Element * > &boundary) const
Definition: nurbs.cpp:2585
NURBSPatchMap(const NURBSExtension *ext)
Definition: nurbs.hpp:460
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:201
int GetNKV() const
Definition: nurbs.hpp:350
int GetNKS() const
Definition: nurbs.hpp:49
int operator()(const int i) const
Definition: nurbs.hpp:601
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:190
Array< int > & GetSlave()
Definition: nurbs.hpp:327
void Generate3DBdrElementDofTable()
Definition: nurbs.cpp:2835
int GetNE() const
Returns number of elements.
Definition: mesh.hpp:846
int NumOfControlPoints
Definition: nurbs.hpp:37
Array< int > e_spaceOffsets
Definition: nurbs.hpp:207
void GenerateBdrElementDofTable()
Definition: nurbs.cpp:2773
const KnotVector * GetKnotVector(int i) const
Definition: nurbs.hpp:363
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:188
void GenerateActiveVertices()
Definition: nurbs.cpp:2001
Array< bool > activeElem
Definition: nurbs.hpp:182
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:359
friend NURBSPatch * Revolve3D(NURBSPatch &patch, double n[], double ang, int times)
Definition: nurbs.cpp:1260
int GetNDof() const
Definition: nurbs.hpp:360
int Dimension() const
Definition: nurbs.hpp:340
Array< int > el_to_patch
Definition: nurbs.hpp:213
Array< int > & GetMaster()
Definition: nurbs.hpp:325
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:206
void UniformRefinement(Vector &newknots) const
Definition: nurbs.cpp:87
const Array< int > & GetSlave() const
Definition: nurbs.hpp:326
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:195
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:153
void Get2DPatchNets(const Vector &Nodes, int vdim)
Definition: nurbs.cpp:3180
Vector & GetWeights()
Definition: nurbs.hpp:382
bool HavePatches() const
Definition: nurbs.hpp:369
virtual ~NURBSExtension()
Destroy a NURBSExtension.
Definition: nurbs.cpp:1677
void GenerateOffsets()
Definition: nurbs.cpp:2318
Array< int > activeDof
Definition: nurbs.hpp:184
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:135
const Array< int > & GetOrders() const
Read-only access to the orders of all knot vectors.
Definition: nurbs.hpp:345
void UniformRefinement()
Definition: nurbs.cpp:671
Array< int > f_spaceOffsets
Definition: nurbs.hpp:208
int GetNBP() const
Definition: nurbs.hpp:342
Array< bool > activeBdrElem
Definition: nurbs.hpp:183
int GetNC() const
Definition: nurbs.hpp:134
void LoadSolution(std::istream &input, GridFunction &sol) const
Definition: nurbs.cpp:3014
int Dimension() const
Definition: mesh.hpp:911
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:348
KnotVector * GetKV(int i)
Definition: nurbs.hpp:136
Array< int > slave
Definition: nurbs.hpp:197
int GetGNV() const
Definition: nurbs.hpp:352
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:241
Array< int > v_meshOffsets
Definition: nurbs.hpp:200
int GetNP() const
Definition: nurbs.hpp:341
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:203
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:371
Table * GetBdrElementDofTable()
Definition: nurbs.hpp:372
Array2D< int > el_to_IJK
Definition: nurbs.hpp:215
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:355
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:209
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:483
void GenerateElementDofTable()
Definition: nurbs.cpp:2623
Array< NURBSPatch * > patches
Definition: nurbs.hpp:218
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
Array< int > mOrders
Definition: nurbs.hpp:173
RefCoord t[3]
Vector data type.
Definition: vector.hpp:60
int operator[](const int i) const
Definition: nurbs.hpp:473
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:357
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:356
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:418
int GetOrder() const
Definition: nurbs.hpp:51
Array< int > bel_to_patch
Definition: nurbs.hpp:214
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