From 2a81f927ca1c21214e58f481aa3c5cfa5594ca69 Mon Sep 17 00:00:00 2001
From: Amaury Johnan <amjohnen@gmail.com>
Date: Fri, 23 Sep 2016 14:18:22 +0000
Subject: [PATCH] When searching reference point from physical point with
Newton-Raphson: start from the reference coordinate of nearer element node.
---
Geo/MElement.cpp | 25 +++++++++++++++++++++++++
Numeric/nodalBasis.h | 3 +++
2 files changed, 28 insertions(+)
diff --git a/Geo/MElement.cpp b/Geo/MElement.cpp
index 24cb8009a1..71af67d7f9 100644
--- a/Geo/MElement.cpp
+++ b/Geo/MElement.cpp
@@ -846,6 +846,31 @@ void MElement::xyz2uvw(double xyz[3], double uvw[3]) const
// routines are implemented for simplices, where the basis functions
// are linear)
uvw[0] = uvw[1] = uvw[2] = 0.;
+
+ // For high order elements, start from the nearer point
+ if (getPolynomialOrder() > 2) {
+ int numNearer = 0;
+ const MVertex *v = getShapeFunctionNode(0);
+ double distNearer = (v->x()-xyz[0])*(v->x()-xyz[0]) +
+ (v->y()-xyz[1])*(v->y()-xyz[1]) +
+ (v->z()-xyz[2])*(v->z()-xyz[2]);
+ for (int i = 1; i < getNumShapeFunctions(); i++) {
+ const MVertex *v = getShapeFunctionNode(i);
+ double dist = (v->x()-xyz[0])*(v->x()-xyz[0]) +
+ (v->y()-xyz[1])*(v->y()-xyz[1]) +
+ (v->z()-xyz[2])*(v->z()-xyz[2]);
+ if (dist < distNearer) {
+ numNearer = i;
+ distNearer = dist;
+ }
+ }
+ const nodalBasis *nb = getFunctionSpace();
+ fullMatrix<double> refpnts = nb->getReferenceNodes();
+ uvw[0] = refpnts(numNearer, 0);
+ uvw[1] = refpnts(numNearer, 1);
+ uvw[2] = refpnts(numNearer, 2);
+ }
+
int iter = 1, maxiter = 20;
double error = 1., tol = 1.e-6;
diff --git a/Numeric/nodalBasis.h b/Numeric/nodalBasis.h
index ec61eff8cb..472aed9c16 100644
--- a/Numeric/nodalBasis.h
+++ b/Numeric/nodalBasis.h
@@ -23,6 +23,9 @@ class nodalBasis {
void getReferenceNodes(fullMatrix<double> &nodes) const {
nodes = points;
}
+ const fullMatrix<double>& getReferenceNodes() const {
+ return points;
+ }
void getReferenceNodesForBezier(fullMatrix<double> &nodes) const;
// Basis functions & gradients evaluation
--
GitLab