Skip to content

Commit

Permalink
Add evaluation of single basis function
Browse files Browse the repository at this point in the history 
Add paraview write for boxes represented by matrices
Add nestedness check of knot vector
Add extra constructor for `gsTensorDomainBoundaryIterator`
FIX ID of `gsTensorDomainIterator`
  • Loading branch information
@hverhelst
hverhelst committed May 16, 2024
1 parent 4d2819a commit 674de66
Show file tree
Hide file tree
Showing 9 changed files with 323 additions and 14 deletions.
23 changes: 23 additions & 0 deletions src/gsCore/gsBasis.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -405,6 +405,24 @@ class gsBasis : public gsFunctionSet<T>

/// @}

/**
@brief Computes function data for basis function \a i
This function evaluates the function \a i and its derivatives at
the points \a in and writes them in the corresponding fields of \a out.
Which field to write (and what to compute) is controlled
by the \a out.flags (see also gsFuncData).
The input points \a in are expected to be compatible with the
implementation/representation of the function, i.e. they should
be points inside the domain of definitition of the function
@param[in] i
@param[in] in
@param[out] out
*/
virtual void computeSingle(index_t i, const gsMatrix<T> & in, gsFuncData<T> & out) const;

inline short_t dim() const {return this->domainDim();}

/*
Expand Down Expand Up @@ -695,6 +713,11 @@ class gsBasis : public gsFunctionSet<T>
virtual void evalAllDers_into(const gsMatrix<T> & u, int n,
std::vector<gsMatrix<T> >& result) const;

/// @brief Evaluate the basis function \a i and its derivatives up
/// to order \a n at points \a u into \a result.
virtual void evalAllDersSingle_into(index_t i, const gsMatrix<T> & u, int n,
std::vector<gsMatrix<T> >& result) const;

/// @brief Evaluate the basis function \a i and its derivatives up
/// to order \a n at points \a u into \a result.
virtual void evalAllDersSingle_into(index_t i, const gsMatrix<T> & u,
Expand Down
54 changes: 54 additions & 0 deletions src/gsCore/gsBasis.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -14,6 +14,7 @@
#pragma once

#include <gsCore/gsBasisFun.h>
#include <gsCore/gsFuncData.h>
#include <gsCore/gsDomainIterator.h>
#include <gsCore/gsBoundary.h>
#include <gsCore/gsGeometry.h>
Expand Down Expand Up @@ -174,6 +175,33 @@ void gsBasis<T>::linearCombination_into(const gsMatrix<T> & coefs,
}
}

template <typename T>
void gsBasis<T>::computeSingle(index_t i, const gsMatrix<T> & in,
gsFuncData<T> & out ) const
{
const unsigned flags = out.flags;

out.dim = this->dimensions();

const int md = out.maxDeriv();
if (md != -1)
evalAllDersSingle_into(i, in, md, out.values);

if (flags & NEED_ACTIVE)
out.actives.setConstant(i);

// if ( flags & NEED_DIV )
// convertValue<T>::derivToDiv(out.values[1], out.divs, info());
// if ( flags & NEED_CURL )
// convertValue<T>::derivToCurl(out.values[1], out.curls, info());
if (flags & NEED_LAPLACIAN)
{
const index_t dsz = out.deriv2Size();
out.laplacians.resize(1, in.cols());
out.laplacians.row(0) = out.values[2].middleRows(dsz*0,out.dim.first).colwise().sum();
}
}


template<class T>
inline gsMatrix<T> gsBasis<T>::laplacian(const gsMatrix<T> & u ) const
Expand Down Expand Up @@ -488,6 +516,32 @@ void gsBasis<T>::evalAllDers_into(const gsMatrix<T> & u, int n,
}
}

template<class T>
void gsBasis<T>::evalAllDersSingle_into(index_t i, const gsMatrix<T> & u, int n,
std::vector<gsMatrix<T> >& result) const
{
result.resize(n+1);

switch(n)
{
case 0:
evalSingle_into(i,u, result[0]);
break;
case 1:
evalSingle_into (i,u, result[0]);
derivSingle_into(i,u, result[1]);
break;
case 2:
evalSingle_into (i,u, result[0]);
derivSingle_into (i,u, result[1]);
deriv2Single_into(i,u, result[2]);
break;
default:
GISMO_ERROR("evalAllDers implemented for order up to 2<"<<n<< " for "<<*this);
break;
}
}

template<class T>
void gsBasis<T>::evalAllDersSingle_into(index_t, const gsMatrix<T> &,
int, gsMatrix<T>&) const
Expand Down
51 changes: 47 additions & 4 deletions src/gsIO/gsWriteParaview.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -215,14 +215,39 @@ template<class T>
void gsWriteParaview(gsBasis<T> const& basis, std::string const & fn,
unsigned npts =NS, bool mesh = false);

/// \brief Export gsHBox to paraview files
/// \brief Export Basis functions to paraview files
///
/// \param basis a basis object
/// \param fn filename where paraview file is written
/// \param npts number of points used for sampling each curve
/// \param mesh if true, the parameter mesh is plotted as well
template<class T>
void gsWriteParaview(gsHBox<2,T> & box, std::string const & fn);
void gsWriteParaview(gsBasis<T> const& basis,
std::vector<index_t> const & indices,
std::string const & fn,
unsigned npts =NS, bool mesh = false);

/// \brief Export an element \a box to paraview files
///
/// \param box an element represented by a bounding box
/// \param fn filename where paraview file is written
/// \param value to write
template<class T>
void gsWriteParaview(const gsMatrix<T> & box, std::string const & fn, T value = 0.0);

/// \brief Export an element \a box to paraview files
///
/// \param box an element represented by a bounding box
/// \param fn filename where paraview file is written
template<class T>
void gsWriteParaview(const gsMatrix<T> & box, const gsVector<T> & values, std::string const & fn);

/// \brief Export gsHBox to paraview files
///
/// \param box a gsHBox
/// \param fn filename where paraview file is written
template<class T>
void gsWriteParaview(const gsHBox<2,T> & box, std::string const & fn);

/// \brief Export gsHBox to paraview files
///
Expand All @@ -231,7 +256,7 @@ void gsWriteParaview(gsHBox<2,T> & box, std::string const & fn);
/// \param npts number of points used for sampling each curve
/// \param mesh if true, the parameter mesh is plotted as well
template<class T>
void gsWriteParaview(gsHBoxContainer<2,T> & box, std::string const & fn);
void gsWriteParaview(const gsHBoxContainer<2,T> & box, std::string const & fn);


/// \brief Export 2D Point set to Paraview file
Expand All @@ -256,6 +281,20 @@ void gsWriteParaviewPoints(gsMatrix<T> const& X,
gsMatrix<T> const& Z,
std::string const & fn);

/// \brief Export 3D Point set to Paraview file
///
/// \param X 1 times n matrix of values for x direction
/// \param Y 1 times n matrix of values for y direction
/// \param Z 1 times n matrix of values for z-direction
/// \param V 1 times n matrix of values for values
/// \param fn filename where paraview file is written
template<class T>
void gsWriteParaviewPoints(gsMatrix<T> const& X,
gsMatrix<T> const& Y,
gsMatrix<T> const& Z,
gsMatrix<T> const& V,
std::string const & fn);

/// \brief Export Point set to Paraview file
///
/// \param points matrix that contain 2D or 3D points, points are columns
Expand Down Expand Up @@ -404,9 +443,13 @@ template<class T>
void writeSingleCompMesh(const gsBasis<T> & basis, const gsGeometry<T> & Geo,
std::string const & fn, unsigned resolution = 8);

/// Export an element \a box
template<class T>
void writeSingleBox(const gsMatrix<T> & box, std::string const & fn, T value);

/// Export a gsHBox
template<class T>
void writeSingleHBox(gsHBox<2,T> & box, std::string const & fn);
void writeSingleHBox(const gsHBox<2,T> & box, std::string const & fn);

/// Export a control net
template<class T>
Expand Down
101 changes: 95 additions & 6 deletions src/gsIO/gsWriteParaview.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -995,7 +995,6 @@ void gsWriteParaview_basisFnct(int i, gsBasis<T> const& basis, std::string const
gsMatrix<T> pts = gsPointGrid(a,b,np) ;

gsMatrix<T> eval_geo = basis.evalSingle ( i, pts ) ;

if ( 3 - d > 0 )
{
np.conservativeResize(3);
Expand Down Expand Up @@ -1175,9 +1174,99 @@ void gsWriteParaview(gsBasis<T> const& basis, std::string const & fn,
collection.save();
}

/// Export Basis functions
template<class T>
void gsWriteParaview(gsBasis<T> const& basis,
const std::vector<index_t> & indices,
std::string const & fn,
unsigned npts, bool mesh)
{
gsParaviewCollection collection(fn);

for (typename std::vector<index_t>::const_iterator idx = indices.cbegin();
idx != indices.cend();
idx++)
{
std::string fileName = fn + util::to_string(*idx);
std::string fileName_nopath = gsFileManager::getFilename(fileName);
gsWriteParaview_basisFnct<T>(*idx, basis, fileName, npts ) ;
collection.addPart(fileName_nopath + ".vts");
}

if ( mesh )
{
std::string fileName = fn + "_mesh";
std::string fileName_nopath = gsFileManager::getFilename(fileName);
writeSingleBasisMesh(basis, fileName);
//collection.addPart(fileName, ".vtp");
collection.addPart(fileName_nopath + ".vtu");
}

collection.save();
}

/// Writes a single \ref gsHBox \a box to a file with name \a fn
template<class T>
void writeSingleBox(const gsMatrix<T> & box, std::string const & fn, T value)
{
gsMatrix<T> points;
gsVector<unsigned> np(box.rows());
np.setConstant(2);
points = gsPointGrid<T>(box.col(0),box.col(1),np);
// The following is needed since gsPointGrid uses gsVector<unsigned> and gsWriteParaviewTPgrid uses gsVector<index_t>...
gsVector<index_t> np2(box.rows());
np2.setConstant(2);


gsMatrix<T> values(1,np.prod());
values.setConstant(value);
gsWriteParaviewTPgrid(points,values,np2,fn);
}

/// Writes \a boxes to a file with name \a fn
template<class T>
void gsWriteParaview(const gsMatrix<T> & boxes, std::string const & fn, T value)
{
gsParaviewCollection collection(fn);

std::string fileName;
for (index_t k=0; k!=boxes.cols()/2; k++)
{
gsMatrix<> tmpbox = boxes.middleCols(2*k,2);
fileName = fn + "_" + util::to_string(k);
writeSingleBox(tmpbox,fileName,value);
fileName = gsFileManager::getFilename(fileName);
collection.addPart(fileName + ".vts");
}

// Write out the collection file
collection.save();
}

/// Writes \a boxes to a file with name \a fn
template<class T>
void gsWriteParaview(const gsMatrix<T> & boxes, const gsVector<T> & values, std::string const & fn)
{
gsParaviewCollection collection(fn);

std::string fileName;
GISMO_ASSERT(boxes.cols()/2==values.rows(),"Number of boxes and values does not match!");
for (index_t k=0; k!=boxes.cols()/2; k++)
{
gsMatrix<> tmpbox = boxes.middleCols(2*k,2);
fileName = fn + "_" + util::to_string(k);
writeSingleBox(tmpbox,fileName,values[k]);
fileName = gsFileManager::getFilename(fileName);
collection.addPart(fileName + ".vts");
}

// Write out the collection file
collection.save();
}

/// Writes a single \ref gsHBox \a box to a file with name \a fn
template<class T>
void writeSingleHBox(gsHBox<2,T> & box, std::string const & fn)
void writeSingleHBox(const gsHBox<2,T> & box, std::string const & fn)
{
gsMatrix<T> points, values(3,4),corners(2,2);
gsVector<index_t> np(2);
Expand All @@ -1192,7 +1281,7 @@ void writeSingleHBox(gsHBox<2,T> & box, std::string const & fn)

/// Writes a single \ref gsHBox \a box to a file with name \a fn
template<class T>
void gsWriteParaview(gsHBox<2,T> & box, std::string const & fn)
void gsWriteParaview(const gsHBox<2,T> & box, std::string const & fn)
{
gsParaviewCollection collection(fn);

Expand All @@ -1205,14 +1294,14 @@ void gsWriteParaview(gsHBox<2,T> & box, std::string const & fn)

/// Writes a container of \ref gsHBox , i.e. a \gsHBoxContainer \a boxes, to a file with name \a fn
template<class T>
void gsWriteParaview(gsHBoxContainer<2,T> & boxes, std::string const & fn)
void gsWriteParaview(const gsHBoxContainer<2,T> & boxes, std::string const & fn)
{
gsParaviewCollection collection(fn);

index_t i=0;
std::string fileName;
for (typename gsHBoxContainer<2,T>::HIterator Hit = boxes.begin(); Hit!=boxes.end(); Hit++)
for (typename gsHBoxContainer<2,T>::Iterator Cit = Hit->begin(); Cit!=Hit->end(); Cit++, i++)
for (typename gsHBoxContainer<2,T>::cHIterator Hit = boxes.cbegin(); Hit!=boxes.cend(); Hit++)
for (typename gsHBoxContainer<2,T>::cIterator Cit = Hit->cbegin(); Cit!=Hit->cend(); Cit++, i++)
{
fileName = fn + util::to_string(i);
writeSingleHBox<T>(*Cit,fileName);
Expand Down
21 changes: 18 additions & 3 deletions src/gsIO/gsWriteParaview_.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -57,10 +57,22 @@ void gsWriteParaview(gsBasis<T> const& basis, std::string const & fn,
unsigned npts, bool mesh);

TEMPLATE_INST
void gsWriteParaview(gsHBox<2,T> & hbox, std::string const & fn);
void gsWriteParaview(gsBasis<T> const& basis,
std::vector<index_t> const & indices,
std::string const & fn,
unsigned npts, bool mesh);

TEMPLATE_INST
void gsWriteParaview(const gsMatrix<T> & box, std::string const & fn, T value);

TEMPLATE_INST
void gsWriteParaview(const gsMatrix<T> & box, const gsVector<T> & values, std::string const & fn);

TEMPLATE_INST
void gsWriteParaview(const gsHBox<2,T> & hbox, std::string const & fn);

TEMPLATE_INST
void gsWriteParaview(gsHBoxContainer<2,T> & hbox, std::string const & fn);
void gsWriteParaview(const gsHBoxContainer<2,T> & hbox, std::string const & fn);

TEMPLATE_INST
void gsWriteParaviewPoints(gsMatrix<T> const& X, gsMatrix<T> const& Y, std::string const & fn);
Expand Down Expand Up @@ -149,7 +161,10 @@ void writeSingleCompMesh(const gsBasis<T> & basis, const gsGeometry<T> & Geo,
std::string const & fn, unsigned resolution);

TEMPLATE_INST
void writeSingleHBox(gsHBox<2,T> & box, std::string const & fn);
void writeSingleBox(const gsMatrix<T> & box, std::string const & fn, T value);

TEMPLATE_INST
void writeSingleHBox(const gsHBox<2,T> & box, std::string const & fn);

TEMPLATE_INST
void writeSingleControlNet(const gsGeometry<T> & Geo, std::string const & fn);
Expand Down
Loading

0 comments on commit 674de66

Please sign in to comment.