Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Adoption of new domain iterators in G+Smo #26

Merged
merged 3 commits into from
Feb 26, 2025
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
20 changes: 11 additions & 9 deletions src/gsC1SurfBasisEdge.h
Original file line number Diff line number Diff line change
Expand Up @@ -263,22 +263,24 @@ namespace gismo
const gsGeometry<T> & patch = m_mp.patch(0);

// Initialize domain element iterator
typename gsBasis<T>::domainIter domIt = m_geo.basis(0).makeDomainIterator(boundary::none);

#ifdef _OPENMP
for ( domIt->next(tid); domIt->good(); domIt->next(nt) )
#else
for (; domIt->good(); domIt->next() )
#endif
typename gsBasis<T>::domainIter domIt = m_geo.basis(0).domain()->beginBdr(boundary::none);
typename gsBasis<T>::domainIter domItEnd = m_geo.basis(0).domain()->endBdr(boundary::none);

# ifdef _OPENMP
domIt += tid;
for ( ; domIt<domItEnd; domIt+=nt )
# else
for (; domIt<domItEnd; ++domIt )
# endif
{
// Map the Quadrature rule to the element
quRule.mapTo( domIt->lowerCorner(), domIt->upperCorner(), quNodes, quWeights );
quRule.mapTo( domIt.lowerCorner(), domIt.upperCorner(), quNodes, quWeights );

// Perform required evaluations on the quadrature nodes
visitor_.evaluate(bf_index, typeBf, basis_g1, basis_geo, basis_plus, basis_minus, patch, quNodes, m_uv, m_gD, m_isBoundary);

// Assemble on element
visitor_.assemble(*domIt, quWeights);
visitor_.assemble(domIt, quWeights);

// Push to global matrix and right-hand side vector
#pragma omp critical(localToGlobal)
Expand Down
20 changes: 11 additions & 9 deletions src/gsC1SurfBasisVertex.h
Original file line number Diff line number Diff line change
Expand Up @@ -230,22 +230,24 @@ namespace gismo
const gsGeometry<T> & patch = m_mp.patch(0);

// Initialize domain element iterator
typename gsBasis<T>::domainIter domIt = m_geo.basis(0).makeDomainIterator(boundary::none);

#ifdef _OPENMP
for ( domIt->next(tid); domIt->good(); domIt->next(nt) )
#else
for (; domIt->good(); domIt->next() )
#endif
typename gsBasis<T>::domainIter domIt = m_geo.basis(0).domain()->beginBdr(boundary::none);
typename gsBasis<T>::domainIter domItEnd = m_geo.basis(0).domain()->endBdr(boundary::none);

# ifdef _OPENMP
domIt += tid;
for ( ; domIt<domItEnd; domIt+=nt )
# else
for (; domIt<domItEnd; ++domIt )
# endif
{
// Map the Quadrature rule to the element
quRule.mapTo( domIt->lowerCorner(), domIt->upperCorner(), quNodes, quWeights );
quRule.mapTo( domIt.lowerCorner(), domIt.upperCorner(), quNodes, quWeights );
#pragma omp critical(evaluate)
// Perform required evaluations on the quadrature nodes
visitor_.evaluate(basis_g1, basis_geo, m_basis_plus, m_basis_minus, patch, quNodes, m_gD, m_isBoundary, m_Phi);

// Assemble on element
visitor_.assemble(*domIt, quWeights);
visitor_.assemble(domIt, quWeights);

// Push to global matrix and right-hand side vector
#pragma omp critical(localToGlobal)
Expand Down
42 changes: 23 additions & 19 deletions src/gsC1SurfGluingData.h
Original file line number Diff line number Diff line change
Expand Up @@ -176,23 +176,25 @@ class gsC1SurfGluingData : public gsC1SurfGD<T>
// Initialize reference quadrature rule and visitor data
visitor_.initialize(basis,quRule);

// Initialize domain element iterator -- using unknown 0
typename gsBasis<T>::domainIter domIt = basis.makeDomainIterator(boundary::none);

#ifdef _OPENMP
for ( domIt->next(tid); domIt->good(); domIt->next(nt) )
#else
for (; domIt->good(); domIt->next() )
#endif
// Initialize domain element iterator
typename gsBasis<T>::domainIter domIt = basis.domain()->beginBdr(boundary::none);
typename gsBasis<T>::domainIter domItEnd = basis.domain()->endBdr(boundary::none);

# ifdef _OPENMP
domIt += tid;
for ( ; domIt<domItEnd; domIt+=nt )
# else
for (; domIt<domItEnd; ++domIt )
# endif
{
// Map the Quadrature rule to the element
quRule.mapTo( domIt->lowerCorner(), domIt->upperCorner(), quNodes, quWeights );
quRule.mapTo( domIt.lowerCorner(), domIt.upperCorner(), quNodes, quWeights );

// Perform required evaluations on the quadrature nodes
visitor_.evaluate(quNodes, this->m_mp);

// Assemble on element
visitor_.assemble(*domIt, quWeights);
visitor_.assemble(domIt, quWeights);

// Push to global matrix and right-hand side vector
#pragma omp critical(localToGlobal)
Expand Down Expand Up @@ -223,23 +225,25 @@ class gsC1SurfGluingData : public gsC1SurfGD<T>
// Initialize reference quadrature rule and visitor data
visitor_Beta.initialize(basis,quRule);

// Initialize domain element iterator -- using unknown 0
typename gsBasis<T>::domainIter domIt = basis.makeDomainIterator(boundary::none);
// Initialize domain element iterator
typename gsBasis<T>::domainIter domIt = basis.domain()->beginBdr(boundary::none);
typename gsBasis<T>::domainIter domItEnd = basis.domain()->endBdr(boundary::none);

#ifdef _OPENMP
for ( domIt->next(tid); domIt->good(); domIt->next(nt) )
#else
for (; domIt->good(); domIt->next() )
#endif
# ifdef _OPENMP
domIt += tid;
for ( ; domIt<domItEnd; domIt+=nt )
# else
for (; domIt<domItEnd; ++domIt )
# endif
{
// Map the Quadrature rule to the element
quRule.mapTo( domIt->lowerCorner(), domIt->upperCorner(), quNodes, quWeights );
quRule.mapTo( domIt.lowerCorner(), domIt.upperCorner(), quNodes, quWeights );

// Perform required evaluations on the quadrature nodes
visitor_Beta.evaluateBeta(quNodes, this->m_mp, sol);

// Assemble on element
visitor_Beta.assembleBeta(*domIt, quWeights);
visitor_Beta.assembleBeta(domIt, quWeights);

// Push to global matrix and right-hand side vector
#pragma omp critical(localToGlobal)
Expand Down
22 changes: 12 additions & 10 deletions src/gsC1SurfGluingDataAssembler.h
Original file line number Diff line number Diff line change
Expand Up @@ -127,7 +127,7 @@ class gsC1SurfGluingDataAssembler : public gsAssembler<T>
m_system_beta_S.matrix().makeCompressed();

}

template <class T, class bhVisitor>
inline void gsC1SurfGluingDataAssembler<T, bhVisitor>::apply(bhVisitor & visitor,
index_t patchIndex,
Expand Down Expand Up @@ -155,23 +155,25 @@ class gsC1SurfGluingDataAssembler : public gsAssembler<T>

//const gsGeometry<T> & patch = m_geo.patch(patchIndex); // 0 = patchindex

// Initialize domain element iterator -- using unknown 0
typename gsBasis<T>::domainIter domIt = basis.makeDomainIterator(boundary::none);
// Initialize domain element iterator
typename gsBasis<T>::domainIter domIt = basis.domain()->beginBdr(boundary::none);
typename gsBasis<T>::domainIter domItEnd = basis.domain()->endBdr(boundary::none);

#ifdef _OPENMP
for ( domIt->next(tid); domIt->good(); domIt->next(nt) )
#else
for (; domIt->good(); domIt->next() )
#endif
# ifdef _OPENMP
domIt += tid;
for ( ; domIt<domItEnd; domIt+=nt )
# else
for (; domIt<domItEnd; ++domIt )
# endif
{
// Map the Quadrature rule to the element
quRule.mapTo( domIt->lowerCorner(), domIt->upperCorner(), quNodes, quWeights );
quRule.mapTo( domIt.lowerCorner(), domIt.upperCorner(), quNodes, quWeights );

// Perform required evaluations on the quadrature nodes
visitor_.evaluate(basis, quNodes, m_uv, m_mp, m_gamma, m_isBoundary);

// Assemble on element
visitor_.assemble(*domIt, quWeights);
visitor_.assemble(domIt, quWeights);

// Push to global matrix and right-hand side vector
#pragma omp critical(localToGlobal)
Expand Down
4 changes: 2 additions & 2 deletions src/gsC1SurfGluingDataVisitor.h
Original file line number Diff line number Diff line change
Expand Up @@ -263,7 +263,7 @@ class gsC1SurfGluingDataVisitor
localRhsBeta.setZero(numActiveBeta, 1 ); //multiple right-hand sides
}

inline void assemble(gsDomainIterator<T> & element,
inline void assemble(gsDomainIteratorWrapper<T> & element,
const gsVector<T> & quWeights)
{
GISMO_UNUSED(element);
Expand All @@ -285,7 +285,7 @@ class gsC1SurfGluingDataVisitor
}
}

inline void assembleBeta(gsDomainIterator<T> & element,
inline void assembleBeta(gsDomainIteratorWrapper<T> & element,
const gsVector<T> & quWeights)
{
GISMO_UNUSED(element);
Expand Down
2 changes: 1 addition & 1 deletion src/gsC1SurfVisitorBasisEdge.h
Original file line number Diff line number Diff line change
Expand Up @@ -158,7 +158,7 @@ namespace gismo
} // Patch 1
} // evaluate1

inline void assemble(gsDomainIterator<T> & element,
inline void assemble(gsDomainIteratorWrapper<T> & element,
const gsVector<T> & quWeights)
{
GISMO_UNUSED(element);
Expand Down
2 changes: 1 addition & 1 deletion src/gsC1SurfVisitorBasisVertex.h
Original file line number Diff line number Diff line change
Expand Up @@ -419,7 +419,7 @@ namespace gismo

} // evaluate

inline void assemble(gsDomainIterator<T> & element,
inline void assemble(gsDomainIteratorWrapper<T> & element,
const gsVector<T> & quWeights)
{
GISMO_UNUSED(element);
Expand Down
43 changes: 25 additions & 18 deletions src/gsContainerBasis.h
Original file line number Diff line number Diff line change
Expand Up @@ -52,40 +52,40 @@ namespace gismo

public:

GISMO_CLONE_FUNCTION(gsContainerBasis)
GISMO_OVERRIDE_CLONE_FUNCTION(gsContainerBasis)

short_t domainDim() const
short_t domainDim() const override
{
return d;
}

void connectivity(const gsMatrix<T> & nodes,
gsMesh<T> & mesh) const
gsMesh<T> & mesh) const override
{
GISMO_UNUSED(nodes); GISMO_UNUSED(mesh);
GISMO_NO_IMPLEMENTATION;
}

memory::unique_ptr<gsGeometry<T> > makeGeometry( gsMatrix<T> coefs ) const
memory::unique_ptr<gsGeometry<T> > makeGeometry( gsMatrix<T> coefs ) const override
{
GISMO_UNUSED(coefs);
GISMO_NO_IMPLEMENTATION;
}

std::ostream &print(std::ostream &os) const
std::ostream &print(std::ostream &os) const override
{
GISMO_UNUSED(os);
GISMO_NO_IMPLEMENTATION;
}

void uniformRefine(int numKnots = 1, int mul=1, int dir=-1)
void uniformRefine(int numKnots = 1, int mul=1, int dir=-1) override
{
for (size_t i=0; i< basisContainer.size(); ++i)
basisContainer[i].uniformRefine(numKnots,mul,dir);
}

// Returm max degree of all the spaces, otherwise i =
short_t degree(short_t dir) const
short_t degree(short_t dir) const override
{
short_t deg = 0;
for (size_t i=0; i< basisContainer.size(); ++i)
Expand All @@ -95,14 +95,14 @@ namespace gismo
return deg;
}

index_t size() const {
index_t size() const override {
index_t sz = 0;
for (size_t i=0; i< basisContainer.size(); ++i)
sz += basisContainer[i].size();
return sz;
}

void reverse()
void reverse() override
{
for (size_t i=0; i< basisContainer.size(); ++i)
{
Expand All @@ -128,14 +128,21 @@ namespace gismo
}
}

index_t nPieces() const {return basisContainer.size();}
index_t nPieces() const override {return basisContainer.size();}

virtual memory::shared_ptr<gsDomain<T> > domain() const override
{
return basisContainer[0].domain();
}

GISMO_DEPRECATED
typename gsBasis<T>::domainIter makeDomainIterator(const boxSide & side) const override
{
// Using the inner basis for iterating
return basisContainer[0].makeDomainIterator(side);
}

GISMO_DEPRECATED
typename gsBasis<T>::domainIter makeDomainIterator() const override
{
// Using the inner basis for iterating
Expand Down Expand Up @@ -166,7 +173,7 @@ namespace gismo
return basisContainer[0].support(i);
}

gsBasis<T>* boundaryBasis_impl(boxSide const & s) const
gsBasis<T>* boundaryBasis_impl(boxSide const & s) const override
{
/*
if (basisContainer.size() != 1)
Expand All @@ -187,7 +194,7 @@ namespace gismo
return bBasis;
}

void active_into(const gsMatrix<T> & u, gsMatrix<index_t> & result) const
void active_into(const gsMatrix<T> & u, gsMatrix<index_t> & result) const override
{
GISMO_ASSERT(u.rows() == d, "Dimension of the points in active_into is wrong");
//GISMO_ASSERT(u.cols() == 1, "Active_into is wrong");
Expand Down Expand Up @@ -215,7 +222,7 @@ namespace gismo
}
}

void eval_into(const gsMatrix<T> & u, gsMatrix<T> & result) const
void eval_into(const gsMatrix<T> & u, gsMatrix<T> & result) const override
{
result.resize(0, u.cols());
for (size_t i=0; i< basisContainer.size(); ++i)
Expand All @@ -228,7 +235,7 @@ namespace gismo
}
}

void deriv_into(const gsMatrix<T> & u, gsMatrix<T> & result) const
void deriv_into(const gsMatrix<T> & u, gsMatrix<T> & result) const override
{
result.resize(0, u.cols());
for (size_t i=0; i< basisContainer.size(); ++i)
Expand All @@ -241,7 +248,7 @@ namespace gismo
}
}

void deriv2_into(const gsMatrix<T> & u, gsMatrix<T> & result) const
void deriv2_into(const gsMatrix<T> & u, gsMatrix<T> & result) const override
{
result.resize(0, u.cols());
for (size_t i=0; i< basisContainer.size(); ++i)
Expand Down Expand Up @@ -284,7 +291,7 @@ namespace gismo
}


gsMatrix<index_t> boundaryOffset(const boxSide & bside, index_t offset) const
gsMatrix<index_t> boundaryOffset(const boxSide & bside, index_t offset) const override
{
gsMatrix<index_t> result;

Expand Down Expand Up @@ -345,7 +352,7 @@ namespace gismo
return result;
}

index_t functionAtCorner(const boxCorner & corner) const
index_t functionAtCorner(const boxCorner & corner) const override
{
if (basisContainer.size() != 1)
{
Expand All @@ -367,7 +374,7 @@ namespace gismo
// basisContainer:
void setBasis(index_t row, gsTensorBSplineBasis<d,T> basis) { basisContainer[row] = basis; }

const gsBasis<T> & piece(const index_t k) const { return basisContainer[k]; }
const gsBasis<T> & piece(const index_t k) const override { return basisContainer[k]; }
// basisContainer END

//std::vector<gsTensorBSplineBasis<d,T>> & getBasisContainer() { return basisContainer; }
Expand Down
Loading