Documentation ¶
Overview ¶
Package gm (geometry and meshes) implements auxiliary functions for handling geometry and mesh structures
Index ¶
- Variables
- func DistPointBoxN(p *PointN, b *BoxN) (dist float64)
- func DistPointLine(p, a, b *Point, tol float64, verbose bool) float64
- func DistPointPoint(a, b *Point) float64
- func DistPointPointN(p *PointN, q *PointN) (dist float64)
- func IsPointIn(p *Point, cMin, cMax []float64, tol float64) bool
- func IsPointInLine(p, a, b *Point, zero, told, tolin float64) bool
- func MetisAdjacency(edges [][2]int, shares map[int][]int) (xadj, adjncy []int32)
- func MetisPartition(edges [][2]int, npart int, recursive bool) (shares map[int][]int, objval int32, parts []int32)
- func MetisPartitionLowLevel(npart, nvert int, xadj, adjncy []int32, recursive bool) (objval int32, parts []int32)
- func MetisShares(edges [][2]int) (shares map[int][]int)
- func PointsLims(pp []*Point) (cmin, cmax []float64)
- func VecDot(u, v []float64) float64
- func VecNew(m float64, u []float64) []float64
- func VecNewAdd(α float64, u []float64, β float64, v []float64) []float64
- func VecNorm(u []float64) float64
- type BezierQuad
- type Bin
- type BinEntry
- type Bins
- func (o *Bins) Append(x []float64, id int, extra interface{})
- func (o Bins) CalcIndex(x []float64) int
- func (o *Bins) Clear()
- func (o Bins) FindAlongSegment(xi, xf []float64, tol float64) []int
- func (o Bins) FindBinByIndex(idx int) *Bin
- func (o Bins) FindClosest(x []float64) (idClosest int, sqDistMin float64)
- func (o *Bins) FindClosestAndAppend(nextID *int, x []float64, extra interface{}, radTol float64, ...) (id int, existent bool)
- func (o *Bins) GetLimits(idxBin int) (xmin, xmax []float64)
- func (o *Bins) Init(xmin, xmax []float64, ndiv []int)
- func (o *Bins) Nactive() (nactive int)
- func (o *Bins) Nentries() (nentries int)
- func (o Bins) String() string
- func (o *Bins) Summary() (l string)
- type BoxN
- type Bspline
- func (o *Bspline) CalcBasis(t float64)
- func (o *Bspline) CalcBasisAndDerivs(t float64)
- func (o *Bspline) Elements() (spans [][]int)
- func (o *Bspline) GetBasis(i int) float64
- func (o *Bspline) GetDeriv(i int) float64
- func (o *Bspline) NumBasis() int
- func (o *Bspline) Point(t float64, option int) (C []float64)
- func (o *Bspline) RecursiveBasis(t float64, i int) float64
- func (o *Bspline) SetControl(Q [][]float64)
- func (o *Bspline) SetOrder(p int)
- type Entity
- type Grid
- func (o *Grid) Boundary(tag int) []int
- func (o *Grid) ContraBasis(m, n, p, k int) la.Vector
- func (o *Grid) ContraMatrix(m, n, p int) *la.Matrix
- func (o *Grid) CovarBasis(m, n, p, k int) la.Vector
- func (o *Grid) CovarMatrix(m, n, p int) *la.Matrix
- func (o *Grid) DetCovarMatrix(m, n, p int) float64
- func (o *Grid) Edge(iEdge int) []int
- func (o *Grid) EdgeGivenTag(tag int) []int
- func (o *Grid) Face(iFace int) []int
- func (o *Grid) FaceGivenTag(tag int) []int
- func (o *Grid) GammaS(m, n, p, k, i, j int) float64
- func (o *Grid) IndexItoMNP(I int) (m, n, p int)
- func (o *Grid) IndexMNPtoI(m, n, p int) (I int)
- func (o *Grid) Lcoeff(m, n, p, k int) float64
- func (o *Grid) MapMeshgrid2d(v la.Vector) (V [][]float64)
- func (o *Grid) MapMeshgrid3d(v la.Vector) (V [][][]float64)
- func (o *Grid) Meshgrid2d() (X, Y [][]float64)
- func (o *Grid) Meshgrid3d() (X, Y, Z [][][]float64)
- func (o *Grid) Ndim() int
- func (o *Grid) Node(I int) (x la.Vector)
- func (o *Grid) Npts(idim int) int
- func (o *Grid) RectGenUniform(xmin, xmax []float64, npts []int)
- func (o *Grid) RectSet2d(X, Y []float64)
- func (o *Grid) RectSet2dU(xmin, xmax, R, S []float64)
- func (o *Grid) RectSet3d(X, Y, Z []float64)
- func (o *Grid) RectSet3dU(xmin, xmax, R, S, T []float64)
- func (o *Grid) SetNurbsSolid(nrb *Nurbs, R, S, T []float64)
- func (o *Grid) SetNurbsSurf2d(nrb *Nurbs, R, S []float64)
- func (o *Grid) SetTransfinite2d(trf *Transfinite, R, S []float64)
- func (o *Grid) SetTransfinite3d(trf *Transfinite, R, S, T []float64)
- func (o *Grid) Size() int
- func (o *Grid) U(m, n, p int) la.Vector
- func (o *Grid) Umax(idim int) float64
- func (o *Grid) Umin(idim int) float64
- func (o *Grid) UnitNormal(N la.Vector, tag, I int)
- func (o *Grid) X(m, n, p int) la.Vector
- func (o *Grid) Xlen(idim int) float64
- func (o *Grid) Xmax(idim int) float64
- func (o *Grid) Xmin(idim int) float64
- type Metrics
- type Nurbs
- func (o *Nurbs) CalcBasis(u []float64)
- func (o *Nurbs) CalcBasisAndDerivs(u []float64)
- func (o *Nurbs) CloneCtrlsAlongCurve(iAlong, jAt int) (Qnew [][][][]float64)
- func (o *Nurbs) CloneCtrlsAlongSurface(iAlong, jAlong, kat int) (Qnew [][][][]float64)
- func (o *Nurbs) ElemBryLocalInds() (I [][]int)
- func (o *Nurbs) Elements() (spans [][]int)
- func (o *Nurbs) ExtractSurfaces() (surfs []*Nurbs)
- func (o *Nurbs) GetBasisI(I []int) float64
- func (o *Nurbs) GetBasisL(l int) float64
- func (o *Nurbs) GetDerivI(dRdu []float64, I []int)
- func (o *Nurbs) GetDerivL(dRdu []float64, l int)
- func (o *Nurbs) GetElemNumBasis() (npts int)
- func (o *Nurbs) GetLimitsQ() (xmin, xmax []float64)
- func (o *Nurbs) GetQ(i, j, k int) (x []float64)
- func (o *Nurbs) GetQl(l int) (x []float64)
- func (o *Nurbs) GetU(dir int) []float64
- func (o *Nurbs) Gnd() int
- func (o *Nurbs) IndBasis(span []int) (L []int)
- func (o *Nurbs) IndsAlongCurve(iAlong, iSpan0, jAt int) (L []int)
- func (o *Nurbs) IndsAlongSurface(iAlong, jAlong, iSpan0, jSpan0, kat int) (L []int)
- func (o *Nurbs) Krefine(X [][]float64) (O *Nurbs)
- func (o *Nurbs) KrefineN(ndiv int, hughesEtAlPaper bool) *Nurbs
- func (o *Nurbs) NonZeroSpans(dir int) [][]int
- func (o *Nurbs) NumBasis(dir int) int
- func (o *Nurbs) Ord(dir int) int
- func (o *Nurbs) Point(C, u []float64, ndim int)
- func (o *Nurbs) PointAndDerivs(...)
- func (o *Nurbs) PointAndFirstDerivs(dCdu *la.Matrix, C, u []float64, ndim int)
- func (o *Nurbs) RecursiveBasis(u []float64, l int) (res float64)
- func (o *Nurbs) SetControl(verts [][]float64, ctrls []int)
- func (o *Nurbs) SetQ(i, j, k int, x []float64)
- func (o *Nurbs) SetQl(l int, x []float64)
- func (o *Nurbs) U(dir, idx int) float64
- func (o *Nurbs) Udelta(dir int) (Δu float64)
- func (o *Nurbs) UfromR(dir int, r float64) (u float64)
- type NurbsExchangeData
- type NurbsExchangeDataSet
- type NurbsPatch
- type Octree
- type Point
- type PointExchangeData
- type PointN
- type Segment
- type Transfinite
Constants ¶
This section is empty.
Variables ¶
var FactoryNurbs = facNurbsT{}
FactoryNurbs generates NURBS'
var FactoryTfinite = facTfinite{}
FactoryTfinite generates transfinite mappings
var (
XDELZERO = 1e-10 // minimum distance between coordinates; i.e. xmax[i]-xmin[i] mininum
)
consts
Functions ¶
func DistPointBoxN ¶
DistPointBoxN returns the distance of a point to the box
NOTE: If point p lies outside box b, the distance to the nearest point on b is returned. If p is inside b or on its surface, zero is returned.
func DistPointLine ¶
DistPointLine computes the distance from p to line passing through a -> b
func DistPointPoint ¶
DistPointPoint computes the unsigned distance from a to b
func DistPointPointN ¶
DistPointPointN returns the distance between two PointN
func IsPointInLine ¶
IsPointInLine returns whether p is inside line passing through a and b
func MetisAdjacency ¶ added in v1.2.1
MetisAdjacency returns adjacency list as a compressed storage format for METIS shares is the map returned by MetisShares
func MetisPartition ¶ added in v1.2.1
func MetisPartition(edges [][2]int, npart int, recursive bool) (shares map[int][]int, objval int32, parts []int32)
MetisPartition performs graph partitioning using METIS This function computes the shares, adjacency list, and partition by calling the other 3 functions
func MetisPartitionLowLevel ¶ added in v1.2.1
func MetisPartitionLowLevel(npart, nvert int, xadj, adjncy []int32, recursive bool) (objval int32, parts []int32)
MetisPartitionLowLevel performs graph partitioning using METIS
func MetisShares ¶ added in v1.2.1
MetisShares returns a map of shares owned by vertices i.e. each vertex is shared by a number of edges, so, we return the map vertexId => [edgeIds...] attached to this vertexId
Example:
0 1 0-------1--------2 | | | 2| 3| 4| | | | 3-------4--------5 5 6
Input. edges = {0,1}, {1,2}, {0,3}, {1,4}, {2,5}, {3,4}, {4,5}
Output. shares = 0:(0,2), 1:(0,1,3), 2:(1,4)
3:(2,5), 4:(3,5,6), 5:(4,6)
(notation. vertexID:(firstEdgeID, secondEdgeID))
NOTE: (1) the pairs or triples will have sorted edgeIDs
(2) len(shares) = number_of_vertices
func PointsLims ¶
PointsLims returns the limits of a set of points
Types ¶
type BezierQuad ¶
type BezierQuad struct { // input Q [][]float64 // control points; can be set outside // auxiliary P []float64 // a point on curve }
BezierQuad implements a quadratic Bezier curve
C(t) = (1-t)² Q0 + 2 t (1-t) Q1 + t² Q2 = t² A + 2 t B + Q0 A = Q2 - 2 Q1 + Q0 B = Q1 - Q0
func (*BezierQuad) DistPoint ¶
func (o *BezierQuad) DistPoint(X []float64) float64
DistPoint returns the distance from a point to this Bezier curve It finds the closest projection which is stored in P
func (*BezierQuad) GetControlCoords ¶
func (o *BezierQuad) GetControlCoords() (X, Y, Z []float64)
GetControlCoords returns the coordinates of control points as 1D arrays (e.g. for plotting)
func (*BezierQuad) GetPoints ¶
func (o *BezierQuad) GetPoints(T []float64) (X, Y, Z []float64)
GetPoints returns points along the curve for given parameter values
func (*BezierQuad) Point ¶
func (o *BezierQuad) Point(C []float64, t float64)
Point returns the x-y-z coordinates of a point on Bezier curve
type BinEntry ¶
type BinEntry struct { ID int // object Id X []float64 // entry coordinate (read only) Extra interface{} // any entity attached to this entry }
BinEntry holds data of an entry to bin
type Bins ¶
type Bins struct { Ndim int // space dimension Xmin []float64 // [ndim] left/lower-most point Xmax []float64 // [ndim] right/upper-most point Xdel []float64 // [ndim] the lengths along each direction (whole box) Size []float64 // size of bins Ndiv []int // [ndim] number of divisions along each direction All []*Bin // [nbins] all bins (there will be an extra "ghost" bin along each dimension) // contains filtered or unexported fields }
Bins implements a set of bins holding entries and is used to fast search entries by given coordinates.
func (Bins) CalcIndex ¶
CalcIndex calculates the bin index where the point x is returns -1 if out-of-range
func (Bins) FindAlongSegment ¶
FindAlongSegment gets the ids of entries that lie close to a segment
Note: the initial (xi) and final (xf) points on segment define a bounding box to filter points
func (Bins) FindBinByIndex ¶
FindBinByIndex finds or allocate new bin corresponding to index idx
func (Bins) FindClosest ¶
FindClosest returns the id of the entry whose coordinates are closest to x
idClosest -- the id of the closest entity. return -1 if out-of-range or not found sqDistMin -- the minimum distance (squared) between x and the closest entity in the same bin NOTE: FindClosest does search the whole area. It only locates neighbours in the same bin where the given x is located. So, if there area no entries in the bin containing x, no entry will be found.
func (*Bins) FindClosestAndAppend ¶
func (o *Bins) FindClosestAndAppend(nextID *int, x []float64, extra interface{}, radTol float64, diff func(idOld int, xNew []float64) bool) (id int, existent bool)
FindClosestAndAppend finds closest point and, if not found, append to bins with a new Id
Input: nextId -- is the Id of the next point. Will be incremented if x is a new point to be added. x -- is the point to be added extra -- extra information attached to point radTol -- is the tolerance for the radial distance (i.e. NOT squared) to decide whether a new point will be appended or not. diff -- [optional] a function for further check that the new and an eventual existent points are really different, even after finding that they coincide (within tol) Output: id -- the id attached to x existent -- flag telling if x was found, based on given tolerance
func (*Bins) Init ¶
Init initialise Bins structure
xmin -- [ndim] min/initial coordinates of the whole space (box/cube) xmax -- [ndim] max/final coordinates of the whole space (box/cube) ndiv -- [ndim] number of divisions for xmax-xmin
type BoxN ¶
type BoxN struct { // essential Lo *PointN // lower point Hi *PointN // higher point // auxiliary ID int // an auxiliary identification number }
BoxN implements a box int he N-dimensional space
func NewBoxN ¶
NewBoxN creates a new box with given limiting coordinates
L -- limits [4] or [6]: xmin,xmax, ymin,ymax, {zmin,zmax} optional
type Bspline ¶
type Bspline struct { // essential T []float64 // array of knots: e.g: T = [t_0, t_1, ..., t_{m-1}] // optional Q [][]float64 // control points (has to call SetControl to change this) // contains filtered or unexported fields }
Bspline holds B-spline data
Reference: [1] Piegl L and Tiller W (1995) The NURBS book, Springer, 646p
func NewBspline ¶ added in v1.0.1
NewBspline returns a new B-spline
func (*Bspline) CalcBasis ¶
CalcBasis computes all non-zero basis functions N[i] @ t Note: use GetBasis to get a particular basis function value
func (*Bspline) CalcBasisAndDerivs ¶
CalcBasisAndDerivs computes all non-zero basis functions N[i] and corresponding first order derivatives of basis functions w.r.t t => dR[i]dt @ t Note: use GetBasis to get a particular basis function value
use GetDeriv to get a particular derivative
func (*Bspline) GetBasis ¶
GetBasis returns the basis function N[i] just computed by CalcBasis or CalcBasisAndDerivs
func (*Bspline) GetDeriv ¶
GetDeriv returns the derivative dN[i]dt just computed by CalcBasisAndDerivs
func (*Bspline) NumBasis ¶
NumBasis returns the number (n) of basis functions == number of control points
func (*Bspline) Point ¶
Point returns the x-y-z coordinates of a point on B-spline option = 0 : use CalcBasis
1 : use RecursiveBasis
func (*Bspline) RecursiveBasis ¶
RecursiveBasis computes one particular basis function N[i] recursively (not efficient)
func (*Bspline) SetControl ¶
SetControl sets B-spline control points
type Entity ¶
type Entity interface { }
Entity defines the geometric (or not) entity/element to be stored in the Octree
type Grid ¶ added in v1.0.1
type Grid struct {
// contains filtered or unexported fields
}
Grid implements (2D/3D) rectangular or curvilinear grid. It also holds metrics data related to curvilinear coordinates.
Notation: m,n,p -- indices used for grid points i,j,k -- indices used for dimension (x,y,z) Ex: the covariant vector @ (m,n,p) is: o.mtr[p][n][m].GovG0 the i component of this vector is: o.mtr[p][n][m].GovG0[i] NOTE: (1) the deep3 structure mtr holds data with the outer index corresponding to z i.e. o.mtr[idxZ][idxY][idxX] (2) the reference coordinates of generated rectangular grids are assumed to be -1 ≤ u ≤ +1
func (*Grid) Boundary ¶ added in v1.0.1
Boundary returns a list of indices of nodes on edge or face of boundary
NOTE: will return empty list if tag is not available see EdgeGivenTag() and FaceGivenTag() functions
func (*Grid) ContraBasis ¶ added in v1.1.0
ContraBasis returns the [k] contravariant basis g^{k} = d{u[k]}/d{x} [@ point m,n,p]
func (*Grid) ContraMatrix ¶ added in v1.1.0
ContraMatrix returns contravariant metrics g^ij = g^i ⋅ g^j [@ point m,n,p]
func (*Grid) CovarBasis ¶ added in v1.1.0
CovarBasis returns the [k] covariant basis g_{k} = d{x}/d{u[k]} [@ point m,n,p]
func (*Grid) CovarMatrix ¶ added in v1.1.0
CovarMatrix returns the covariant metrics g_ij = g_i ⋅ g_j [@ point m,n,p]
func (*Grid) DetCovarMatrix ¶ added in v1.1.0
DetCovarMatrix returns the determinant of covariant g matrix = det(CovariantMatrix) [@ point m,n,p]
func (*Grid) Edge ¶ added in v1.0.1
Edge returns the ids of points on edges: [edge0, edge1, edge2, edge3]
3 +-----------+ Considering the x-y axes below, the order of indices follows: | | | | y 0 1 2 3 0| |1 ↑ {xmin_edge, xmax_edge, ymin_edge, ymax_edge} | | │ | | +——→ x +-----------+ 2
func (*Grid) EdgeGivenTag ¶ added in v1.1.0
EdgeGivenTag returns a list of nodes marked with given tag
21 +-----------+ Considering the x-y axes below, the order of tags follows: | | | | y 0 1 2 3 10| |11 ↑ {10, 11, 20, 21} | | │ | | +——→ x +-----------+ 20 NOTE: will return empty list if tag is not available
func (*Grid) Face ¶ added in v1.0.1
Face returns the ids of points on faces: [face0, face1, face2, face3, face4, face5]
+----------------+ Considering the x-y-z axes below, ,'| ,'| the order of indices follows: ,' | ___ ,' | ,' |,'5,' [0] ,' | z 0: xmin_face ,' |~~~ ,' , | ↑ 1: xmax_face +'===============+' ,'| | │ 2: ymin_face | ,'| | | |3| | +——→y 3: ymax_face | |2| | | |,' | ,' 4: zmin_face | |,' +- - - | +- - - -+ x 5: zmax_face | ' ,' | ,' | ,' [1] ___| ,' | ,' ,'4,'| ,' | ,' ~~~ | ,' +----------------+'
func (*Grid) FaceGivenTag ¶ added in v1.1.0
FaceGivenTag returns a list of nodes marked with given tag
+----------------+ Considering the x-y-z axes below, ,'| ,'| the order of tags follows: ,' | ___ ,' | ,' |,301' 100 ,' | z 100: xmin_face ,' |~~~ ,' ,' | ↑ 101: xmax_face +'===============+' ,' | | │ 200: ymin_face | ,'| | | |201 | +——→y 201: ymax_face | |200 | | |,' | ,' 300: zmin_face | | ,' +- - - | +- - - -+ x 301: zmax_face | ,' ,' | ,' | ,'101 ___| ,' | ,' ,300'| ,' | ,' ~~~ | ,' +----------------+' NOTE: will return empty list if tag is not available
func (*Grid) GammaS ¶ added in v1.1.0
GammaS returns the [k][i][j] Christoffel coefficients of second kind [@ point m,n,p]
func (*Grid) IndexItoMNP ¶ added in v1.1.0
IndexItoMNP converts node index I into triplet indices (m,n,p)
2D: I = m + n⋅n0 m = I % n0 n = I / n0 3D: I = m + n⋅n0 + p⋅n0⋅n1 p = I / (n0⋅n1) t = I % (n0⋅n1) (projection @ z=0) m = t % n0 n = t / n0
func (*Grid) IndexMNPtoI ¶ added in v1.1.0
IndexMNPtoI converts node triplet indices (m,n,p) into node index I
2D: I = m + n⋅n0 m = I % n0 n = I / n0 3D: I = m + n⋅n0 + p⋅n0⋅n1 p = I / (n0⋅n1) t = I % (n0⋅n1) (projection @ z=0) m = t % n0 n = t / n0
func (*Grid) Lcoeff ¶ added in v1.1.0
Lcoeff returns the [k] L-coefficients = sum(Γ_ij^k ⋅ g^ij) [@ point m,n,p]
func (*Grid) MapMeshgrid2d ¶ added in v1.1.0
MapMeshgrid2d maps vector V into 2D meshgrid using node indices conversion IndexMNPtoI()
vv[ny][nx] -- mapped values: vv[n][m] ⇐ V[I] (see also Meshgrid2d)
func (*Grid) MapMeshgrid3d ¶ added in v1.1.0
MapMeshgrid3d maps vector V into 3D meshgrid using node indices conversion IndexMNPtoI()
V[nz][ny][nx] -- mapped values: V[p][n][m] ⇐ v[I] (see also Meshgrid2d)
func (*Grid) Meshgrid2d ¶ added in v1.1.0
Meshgrid2d extracts 2D meshgrid
X -- x0[ny][nx] Y -- x1[ny][nx]
func (*Grid) Meshgrid3d ¶ added in v1.1.0
Meshgrid3d extracts 3D meshgrid
X -- x0[nz][ny][nx] Y -- x1[nz][ny][nx] Z -- x2[nz][ny][nx]
func (*Grid) Node ¶ added in v1.0.1
Node returns the physical coordinates of node I. See IndexItoMNP(I) ⇒ (m,n,p).
x -- slice to position vector @ [p][n][m] [may be used to change values]
func (*Grid) RectGenUniform ¶ added in v1.1.0
RectGenUniform generates uniform coordinates of a rectangular grid
xmin -- min x-y-z values, len==ndim: [xmin, ymin, zmin] xmax -- max x-y-z values, len==ndim: [xmax, ymax, zmax] npts -- number of points along each direction len==ndim: [n0, n1, n2] (must be greater than 2) -1 ≤ u ≤ +1 x(u) = xmin + (xmax - xmin) ⋅ (1 + u) / 2 u(x) = -1 + 2⋅(x - xmin) / (xmax - xmin) dx/du = (xmax - xmin) / 2
func (*Grid) RectSet2d ¶ added in v1.1.0
RectSet2d sets rectangular grid with given coordinates
-1 ≤ u ≤ +1 x(u) = xmin + (xmax - xmin) ⋅ (1 + u) / 2 u(x) = -1 + 2⋅(x - xmin) / (xmax - xmin) dx/du = (xmax - xmin) / 2
func (*Grid) RectSet2dU ¶ added in v1.1.0
RectSet2dU sets rectangular grid with given reference coordinates and limits
Input: xmin -- min x-y values: [xmin, ymin] xmax -- max x-y values: [xmax, ymax] R -- reference coordinates along x: -1 ≤ r ≤ +1 S -- reference coordinates along y: -1 ≤ s ≤ +1 -1 ≤ u ≤ +1 x(u) = xmin + (xmax - xmin) ⋅ (1 + u) / 2 u(x) = -1 + 2⋅(x - xmin) / (xmax - xmin) dx/du = (xmax - xmin) / 2
func (*Grid) RectSet3d ¶ added in v1.1.0
RectSet3d sets rectangular grid with given coordinates
-1 ≤ u ≤ +1 x(u) = xmin + (xmax - xmin) ⋅ (1 + u) / 2 u(x) = -1 + 2⋅(x - xmin) / (xmax - xmin) dx/du = (xmax - xmin) / 2
func (*Grid) RectSet3dU ¶ added in v1.1.0
RectSet3dU sets rectangular grid with given reference coordinates and limits
Input: xmin -- min x-y-z values: [xmin, ymin, zmin] xmax -- max x-y-z values: [xmax, ymax, zmax] R -- reference coordinates along x: -1 ≤ r ≤ +1 S -- reference coordinates along y: -1 ≤ s ≤ +1 T -- reference coordinates along z: -1 ≤ t ≤ +1 -1 ≤ u ≤ +1 x(u) = xmin + (xmax - xmin) ⋅ (1 + u) / 2 u(x) = -1 + 2⋅(x - xmin) / (xmax - xmin) dx/du = (xmax - xmin) / 2
func (*Grid) SetNurbsSolid ¶ added in v1.1.0
SetNurbsSolid sets grid with NURBS solid
nrb -- NURBS solid R -- [n0] reference coordinates along r-direction -1 ≤ r ≤ +1 S -- [n1] reference coordinates along s-direction -1 ≤ s ≤ +1 T -- [n2] reference coordinates along s-direction -1 ≤ s ≤ +1
func (*Grid) SetNurbsSurf2d ¶ added in v1.1.0
SetNurbsSurf2d sets grid with NURBS surface in 2D (flat surface)
nrb -- NURBS surface in 2D (flat) R -- [n0] reference coordinates along r-direction -1 ≤ r ≤ +1 S -- [n1] reference coordinates along s-direction -1 ≤ s ≤ +1
func (*Grid) SetTransfinite2d ¶ added in v1.1.0
func (o *Grid) SetTransfinite2d(trf *Transfinite, R, S []float64)
SetTransfinite2d sets grid from 2D transfinite mapping
trf -- 2D transfinite structure R -- [n0] reference coordinates along r-direction -1 ≤ r ≤ +1 S -- [n1] reference coordinates along s-direction -1 ≤ s ≤ +1
func (*Grid) SetTransfinite3d ¶ added in v1.1.0
func (o *Grid) SetTransfinite3d(trf *Transfinite, R, S, T []float64)
SetTransfinite3d sets grid from 3D transfinite mapping
trf -- 2D transfinite structure R -- [n0] reference coordinates along r-direction -1 ≤ r ≤ +1 S -- [n1] reference coordinates along s-direction -1 ≤ s ≤ +1 T -- [n2] reference coordinates along s-direction -1 ≤ t ≤ +1
func (*Grid) UnitNormal ¶ added in v1.1.0
UnitNormal computes the unit normal vector at an edge or face defined by "tag" and at a node specified by index "I".
Input: tag -- tag of edge or face; see EdgeGivenTag() or FaceGivenTag() I -- node index; see IndexMNPtoI() [must be a node on boundary; otherwise, panic] Output: N -- unit normal vector NOTE: this function does not check whether point I is on the selected edge or not
func (*Grid) Xlen ¶ added in v1.1.0
Xlen returns the lengths along each direction (whole box) == Xmax(idim) - Xmin(idim)
type Metrics ¶ added in v1.1.0
type Metrics struct { U la.Vector // reference coordinates {r,s,t} X la.Vector // physical coordinates {x,y,z} CovG0 la.Vector // covariant basis g_0 = d{x}/dr CovG1 la.Vector // covariant basis g_1 = d{x}/ds CovG2 la.Vector // covariant basis g_2 = d{x}/dt CntG0 la.Vector // contravariant basis g_0 = dr/d{x} (gradients) CntG1 la.Vector // contravariant basis g_1 = ds/d{x} (gradients) CntG2 la.Vector // contravariant basis g_2 = dt/d{x} (gradients) CovGmat *la.Matrix // covariant metrics g_ij = g_i ⋅ g_j CntGmat *la.Matrix // contravariant metrics g^ij = g^i ⋅ g^j DetCovGmat float64 // determinant of covariant g matrix = det(CovGmat) Homogeneous bool // homogeneous grid => nil second order derivatives and Christoffel symbols GammaS [][][]float64 // [k][i][j] Christoffel coefficients of second kind (non-homogeneous) L []float64 // [3] L-coefficients = sum(Γ_ij^k ⋅ g^ij) (non-homogeneous) }
Metrics holds data related to a position in a space represented by curvilinear coordinates
func NewMetrics2d ¶ added in v1.1.0
NewMetrics2d allocate new 2D metrics structure
NOTE: the second order derivatives (from ddxdrr) may be nil => homogeneous grid
func NewMetrics3d ¶ added in v1.1.0
func NewMetrics3d(u, x, dxdr, dxds, dxdt, ddxdrr, ddxdss, ddxdtt, ddxdrs, ddxdrt, ddxdst la.Vector) (o *Metrics)
NewMetrics3d allocate new 3D metrics structure
NOTE: the second order derivatives (from ddxdrr) may be nil => homogeneous grid
type Nurbs ¶
type Nurbs struct { Q [][][][]float64 // Qw: weighted control points and weights [n[0]][n[1]][n[2]][4] (Piegl p120) // contains filtered or unexported fields }
Nurbs holds NURBS data
NOTE: (1) Control points must be set after a call to Init (2) Either SetControl must be called or the Q array must be directly specified Reference: [1] Piegl L and Tiller W (1995) The NURBS book, Springer, 646p
func (*Nurbs) CalcBasis ¶
CalcBasis computes all non-zero basis functions R[i][j][k] @ u. Note: use GetBasisI or GetBasisL to get a particular basis function value
func (*Nurbs) CalcBasisAndDerivs ¶
CalcBasisAndDerivs computes all non-zero basis functions R[i][j][k] and corresponding first order derivatives of basis functions w.r.t u => dRdu[i][j][k] @ u Note: use GetBasisI or GetBasisL to get a particular basis function value
use GetDerivI or GetDerivL to get a particular derivative
func (*Nurbs) CloneCtrlsAlongCurve ¶
CloneCtrlsAlongCurve returns a copy of control points @ 2D boundary
func (*Nurbs) CloneCtrlsAlongSurface ¶
CloneCtrlsAlongSurface returns a copy of control points @ 3D boundary
func (*Nurbs) ElemBryLocalInds ¶
ElemBryLocalInds returns the local (element) indices of control points @ boundaries (if element would have all surfaces @ boundaries)
func (*Nurbs) ExtractSurfaces ¶
ExtractSurfaces returns a new NURBS representing a boundary of this NURBS
func (*Nurbs) GetBasisI ¶
GetBasisI returns the basis function R[i][j][k] just computed by CalcBasis or CalcBasisAndDerivs Note: I = [i,j,k]
func (*Nurbs) GetBasisL ¶
GetBasisL returns the basis function R[i][j][k] just computed by CalcBasis or CalcBasisAndDerivs Note: l = i + j * n0 + k * n0 * n1
func (*Nurbs) GetDerivI ¶
GetDerivI returns the derivative of basis function dR[i][j][k]du just computed by CalcBasisAndDerivs Note: I = [i,j,k]
func (*Nurbs) GetDerivL ¶
GetDerivL returns the derivative of basis function dR[i][j][k]du just computed by CalcBasisAndDerivs Note: l = i + j * n0 + k * n0 * n1
func (*Nurbs) GetElemNumBasis ¶
GetElemNumBasis returns the number of control points == basis functions needed for one element
npts := Π_i (p[i] + 1)
func (*Nurbs) GetLimitsQ ¶
GetLimitsQ computes the limits of all coordinates of control points in NURBS
func (*Nurbs) IndBasis ¶
IndBasis returns the indices of basis functions == local indices of control points
func (*Nurbs) IndsAlongCurve ¶
IndsAlongCurve returns the control points indices along curve
func (*Nurbs) IndsAlongSurface ¶
IndsAlongSurface return the control points indices along surface
func (*Nurbs) KrefineN ¶
KrefineN return a new Nurbs with each span divided into ndiv parts = [2, 3, ...]
func (*Nurbs) NonZeroSpans ¶
NonZeroSpans returns the 'elements' along direction dir
func (*Nurbs) Point ¶
Point returns the x-y-z coordinates of a point on curve/surface/solid
Input: u -- [gnd] knot values ndim -- the dimension of the point. E.g. allows drawing curves in 3D Output: C -- [ndim] point coordinates NOTE: Algorithm A4.1 (p124) of [1]
func (*Nurbs) PointAndDerivs ¶ added in v1.1.0
func (o *Nurbs) PointAndDerivs(x, dxdr, dxds, dxdt, ddxdrr, ddxdss, ddxdtt, ddxdrs, ddxdrt, ddxdst, u la.Vector, ndim int)
PointAndDerivs computes position and the first and second order derivatives Using Algorithms A3.2(p93), A3.6(p111), A4.2(p127), and A4.4(p137)
Input: u -- knot values {r,s,t} [gnd] ndim -- the dimension of the point. E.g. allows drawing curves in 3D [ndim=3 if gnd=3] Output: x -- position {x,y,z} (the same as the C varible in [1]) dxdr -- ∂{x}/∂r dxds -- ∂{x}/∂s [may be nil] (solid and surfaces) dxdt -- ∂{x}/∂t [may be nil] (solid) ddxdrr -- ∂²{x}/∂r² [optional] ddxdss -- ∂²{x}/∂s² [optional] [may be nil] (solid and surfaces) ddxdtt -- ∂²{x}/∂t² [optional] [may be nil] (solid) ddxdrs -- ∂²{x}/∂r∂s [optional] [may be nil] (solid and surfaces) ddxdrt -- ∂²{x}/∂r∂t [optional] [may be nil] (solid) ddxdst -- ∂²{x}/∂s∂t [optional] [may be nil] (solid) NOTE: the second order derivatives will be ignored if ddxdrr == nil; otherwise, all 2nd derivs will be computed
func (*Nurbs) PointAndFirstDerivs ¶ added in v1.1.0
PointAndFirstDerivs returns the point and first order derivatives with respect to the knot values u of the x-y-z coordinates of a point on curve/surface/solid
Input: u -- [gnd] knot values ndim -- the dimension of the point. E.g. allows drawing curves in 3D Output: dCdu -- [ndim][gnd] derivatives dC_i/du_j C -- [ndim] point coordinates
func (*Nurbs) RecursiveBasis ¶
RecursiveBasis implements basis functions by means of summing all terms in Bernstein polynomial using recursive Cox-DeBoor formula (very not efficient) Note: l = i + j * n0 + k * n0 * n1
func (*Nurbs) SetControl ¶
SetControl sets control points from list of global vertices
type NurbsExchangeData ¶
type NurbsExchangeData struct { ID int `json:"i"` // id of Nurbs Gnd int `json:"g"` // 1: curve, 2:surface, 3:volume (geometry dimension) Ords []int `json:"o"` // order along each x-y-z direction [gnd] Knots [][]float64 `json:"k"` // knots along each x-y-z direction [gnd][m] Ctrls []int `json:"c"` // global ids of control points }
NurbsExchangeData holds all data required to exchange NURBS; e.g. read/save files
type NurbsExchangeDataSet ¶
type NurbsExchangeDataSet []*NurbsExchangeData
NurbsExchangeDataSet defines a set of nurbs exchange data
type NurbsPatch ¶
type NurbsPatch struct { // input/output data ControlPoints []*PointExchangeData `json:"points"` // input/output control points ExchangeData NurbsExchangeDataSet `json:"patch"` // input/output nurbs data // Nurbs structures Entities []*Nurbs `json:"-"` // pointers to NURBS // auxiliary Bins Bins `json:"-"` // auxiliary structure to locate points }
NurbsPatch defines a patch of many NURBS'
func NewNurbsPatch ¶
func NewNurbsPatch(binsNdiv int, tolerance float64, entities ...*Nurbs) (o *NurbsPatch)
NewNurbsPatch returns new patch of NURBS
tolerance -- tolerance to assume that two control points are the same
func NewNurbsPatchFromFile ¶
func NewNurbsPatchFromFile(filename string, binsNdiv int, tolerance float64) (o *NurbsPatch)
NewNurbsPatchFromFile allocates a NurbsPatch with data from file
func (NurbsPatch) LimitsAndNdim ¶
func (o NurbsPatch) LimitsAndNdim() (xmin, xmax []float64, ndim int)
LimitsAndNdim computes the limits of the patch and max dimension by looping over all Entities
func (*NurbsPatch) ResetFromEntities ¶
func (o *NurbsPatch) ResetFromEntities(binsNdiv int, tolerance float64)
ResetFromEntities will reset all exchange data with information from Entities slice
func (*NurbsPatch) ResetFromExchangeData ¶
func (o *NurbsPatch) ResetFromExchangeData(binsNdiv int, tolerance float64)
ResetFromExchangeData will reset all Entities with information from ExchangeData (and ControlPoints)
func (NurbsPatch) Write ¶
func (o NurbsPatch) Write(dirout, fnkey string)
Write writes ExchangeData to json file
type Octree ¶
type Octree struct { // constants DIM uint32 // dimension PMAX uint32 // roughly how many levels fit in 32 bits QO uint32 // 4 for quadtree, 8 for octree QL uint32 // offset constant to leftmost daughter // contains filtered or unexported fields }
Octree implements a Quad-Tree or an Oct-Tree to assist in fast-searching elements (entities) in the 2D or 3D space
type Point ¶
type Point struct {
X, Y, Z float64
}
Point holds the Cartesian coordinates of a point in 3D space
type PointExchangeData ¶
type PointExchangeData struct { ID int `json:"i"` // id Tag int `json:"t"` // tag X []float64 `json:"x"` // coordinates (size==4) }
PointExchangeData holds data for exchanging control points
type PointN ¶
type PointN struct { // esssential X []float64 // coordinates // optional ID int // some identification number }
PointN implements a point in N-dim space
func NewPointNdim ¶
NewPointNdim creates a new PointN with given dimension (ndim)
func (PointN) AlmostTheSameX ¶
AlmostTheSameX returns true if the X slices of two points have almost the same values, for given tolerance
func (PointN) ExactlyTheSameX ¶
ExactlyTheSameX returns true if the X slices of two points have exactly the same values
type Segment ¶
type Segment struct {
A, B *Point
}
Segment represents a directed segment from A to B
func NewSegment ¶
NewSegment creates a new segment from a to b
type Transfinite ¶ added in v1.0.1
type Transfinite struct {
// contains filtered or unexported fields
}
Transfinite maps a reference square [-1,+1]×[-1,+1] into a curve-bounded quadrilateral
B[3](r(x,y)) _,'\ B[3](r) _,' \ B[1](s(x,y)) ┌───────┐ _,' \ │ │ \ \ B[0](s)│ │B[1](s) ⇒ \ _,' s │ │ B[0](s(x,y)) \ _,' │ └───────┘ y \ _,' B[2](r(x,y)) └──r B[2](r) │ ' └──x +----------------+ ,'| ,'| t or z ,' | ___ ,' | B[0](s,t) ↑ ,' |,'5,' [0],' | B[1](s,t) | ,' |~~~ ,' | B[2](r,t) | +'===============+' ,'| | B[3](r,t) | | ,'| | | |3| | B[4](r,s) | s or y | |2| | | |,' | B[5](r,s) +--------> | |,' +- - - | +- - - -+ ,' | ,' | ,' ,' | ,' [1] ___| ,' r or x | ,' ,'4,'| ,' | ,' ~~~ | ,' +----------------+' NOTE: the "reference" coordinates {r,s,t} ϵ [-1,+1]×[-1,+1]×[-1,+1] are also symbolised as "u"
func NewTransfinite2d ¶ added in v1.1.0
func NewTransfinite2d(B, Bd, Bdd []fun.Vs) (o *Transfinite)
NewTransfinite2d allocates a new structure
Input: B -- [4] boundary functions Bd -- [4] 1st derivative of boundary functions Bdd -- [4 or nil] 2nd derivative of boundary functions [may be nil]
func NewTransfinite3d ¶ added in v1.1.0
NewTransfinite3d allocates a new structure
Input: B -- [6] boundary functions Bd -- [6] 1st derivative of boundary functions Bdd -- [6 or nil] 2nd derivative of boundary functions [may be nil]
func (*Transfinite) Point ¶ added in v1.0.1
func (o *Transfinite) Point(x, u la.Vector)
Point computes "real" position x(r,s,t)
Input: u -- the "reference" coordinates {r,s,t} ϵ [-1,+1]×[-1,+1]×[-1,+1] Output: x -- the "real" coordinates {x,y,z}
func (*Transfinite) PointAndDerivs ¶ added in v1.1.0
func (o *Transfinite) PointAndDerivs(x, dxDr, dxDs, dxDt, ddxDrr, ddxDss, ddxDtt, ddxDrs, ddxDrt, ddxDst, u la.Vector)
PointAndDerivs computes position and the first and second order derivatives
Input: u -- reference coordinates {r,s,t} Output: x -- position {x,y,z} dxDr -- ∂{x}/∂r dxDs -- ∂{x}/∂s dxDt -- ∂{x}/∂t ddxDrr -- ∂²{x}/∂r² [may be nil] ddxDss -- ∂²{x}/∂s² [may be nil] ddxDtt -- ∂²{x}/∂t² [may be nil] ddxDrs -- ∂²{x}/∂r∂s [may be nil] ddxDrt -- ∂²{x}/∂r∂t [may be nil] ddxDst -- ∂²{x}/∂s∂t [may be nil]
Source Files ¶
Directories ¶
Path | Synopsis |
---|---|
Package msh defines mesh data structures and implements interpolation functions for finite element analyses (FEA)
|
Package msh defines mesh data structures and implements interpolation functions for finite element analyses (FEA) |
Package tri wraps Triangle to perform mesh generation a Delaunay triangulation
|
Package tri wraps Triangle to perform mesh generation a Delaunay triangulation |