From: Peter Schaefer <peter.schaefer@tuwien.ac.at>
Date: Wed, 13 Mar 2013 08:37:55 +0000 (+0100)
Subject: [src] QE Bug behoben, Falsche Funktion wurde aufgerufen
X-Git-Url: https://git.leopard-lacewing.eu/?a=commitdiff_plain;h=da8022067e518654c9f55de1102eb285bf0dd75a;p=bacc.git

[src] QE Bug behoben, Falsche Funktion wurde aufgerufen
---

diff --git a/src/slpRectangle.cpp b/src/slpRectangle.cpp
index 5283e4c..962762b 100644
--- a/src/slpRectangle.cpp
+++ b/src/slpRectangle.cpp
@@ -13,7 +13,7 @@
 //#include <iostream>
 #include <math.h>
 //#include <cassert>
-//#include <limits>
+#include <limits>
 //#include <cstdlib>
 
 #include <mex.h>
@@ -51,7 +51,9 @@ double inline max(double x, double y) {
 double inline min(double x, double y) {
 	return x > y ? y : x;
 }
-
+/*
+ * kleinere Seiten nach vorn und dadurch kleinerer Durchmesser nach vorn
+ */
 void inline switch_site(double& b, double& d, double& t, double& v, double& d1,
 		double& d2) {
 	double tmp = 0;
@@ -75,6 +77,7 @@ void inline switch_site(double& b, double& d, double& t, double& v, double& d1,
 
 void inline switch_dim(double& b, double& d, double& t, double& v, double& d1,
 		double& d2) {
+	//
 	if (b * b + t * t > d * d + v * v) {
 		double tmp = 0;
 
@@ -91,6 +94,9 @@ void inline switch_dim(double& b, double& d, double& t, double& v, double& d1,
 	}
 }
 
+/*
+ * gibt den Abstand zum Quadrat zurueck
+ */
 double inline dist2(double b, double d, double t, double v, double d1,
 		double d2, double d3) {
 	double dis1 = 0, dis2 = 0;
@@ -110,11 +116,17 @@ double inline dist2(double b, double d, double t, double v, double d1,
 
 	return dis1 * dis1 + dis2 * dis2 + d3 * d3;
 }
+/*
+ * gibt den Abstand zurueck
+ */
 double inline dist(double b, double d, double t, double v, double d1, double d2,
 		double d3) {
 	return sqrt(dist2(b, d, t, v, d1, d2, d3));
 }
 
+/*
+ * gibt den Abstand in einer Richtung zum Quadrat zurueck
+ */
 double inline dist_s2(double b, double t, double d1) {
 	double dis1 = 0;
 	if (d1 < 0)
@@ -126,6 +138,9 @@ double inline dist_s2(double b, double t, double d1) {
 
 	return dis1 * dis1;
 }
+/*
+ * gibt den Abstand in einer Richtung zurueck
+ */
 double inline dist_s(double b, double t, double d1) {
 	return sqrt(dist_s2(b, t, d1));
 }
@@ -210,11 +225,12 @@ double inline G_AY2X2(double y1, double y2, double x1, double x2, double d3) {
 		//sol -= (x2 - y2) * g_AY(-0.5,x2,y2, sqrt((x1 - y1) * (x1 - y1) + d3 * d3));
 		sol += (x2 - y2) * log(sqrt(hlp) - (x2 - y2));
 
-	if (sol != sol /*|| fabs(sol) == numeric_limits<double>::infinity()*/) {
+	if (sol != sol/* || fabs(sol) == std::numeric_limits<double>::infinity()*/) {
 		//	cout << "G_AY2X2: " << sol << " isn't a Number. (" << y1 << "," << y2
 		//			<< ")(" << x1 << "," << x2 << ")(" << d3 << ")" << endl;
 		//	cout << (x2 - y2) << "__" << sqrt(hlp) << endl;
 		//	cout << (x1 - y1) * (x1 - y1) << " " << hlp << endl;
+		// mexPrintf("G_AY2X2X: nan");
 	}
 	return sol;
 }
@@ -414,8 +430,8 @@ double calcParIntQY1X1(double b, double d, double t, double v, double d1,
 					b * GAUSS_NODES[quad][i].n, 0, d3) * t * b
 			* GAUSS_NODES[quad][i].w * GAUSS_NODES[quad][j].w;
 
-			if (sol
-					!= sol /*|| fabs(sol) == numeric_limits<double>::infinity()*/) {
+			if (sol!= sol /*|| fabs(sol) == numeric_limits<double>::infinity()*/) {
+				//mexPrintf("(%f,%f),(%f,%f),(%f,%f,%f)",b,d,t,v,d1,d2,d3);
 				/*				cout << "->(" << i << "," << j << ")" << endl;
 				 cout << "-|("
 				 << G_AY2X2(t * GAUSS_NODES[quad][j].n + d1, v + d2,
@@ -446,9 +462,9 @@ double calcParIntQY1X1(double b, double d, double t, double v, double d1,
 				 * GAUSS_NODES[quad][j].w << ")" << endl;
 
 				 cout << endl;*/
-				return sol;
+				//mexPrintf("%f\n",sol);
 			}
-
+			return sol;
 		}
 	}
 
@@ -546,7 +562,10 @@ double cParO2(double b, double d, double t, double v, double d1, double d2,
 	//kurze Seite nach vorn
 	switch_site(b, d, t, v, d1, d2);
 
-	if ((t * t + v * v) <= zeta[0] * dist2(b, d, t, v, d1, d2, d3))
+	//kurze Achse nach vorn
+//	switch_dim(b, d, t, v, d1, d2);
+
+	if (zeta[0]*zeta[0] *(t * t + v * v) <  dist2(b, d, t, v, d1, d2, d3))
 		return calcParIntQ(b, d, t, v, d1, d2, d3);
 
 	return calcParIntA(b, d, t, v, d1, d2, d3);
@@ -559,15 +578,15 @@ double cParO3(double b, double d, double t, double v, double d1, double d2,
 	switch_site(b, d, t, v, d1, d2);
 
 	//kurze Achse nach vorn
-	switch_dim(b, d, t, v, d1, d2);
+//	switch_dim(b, d, t, v, d1, d2);
 
-	if ((b * b + d * d) <= zeta[1] * dist2(b, d, t, v, d1, d2, d3))
-		if ((max(b, t) <= zeta[2] * min(d, v))
-				&& (min(dist_s2(b, t, d1), dist_s2(d, v, d2)) > 0))
+	if (zeta[1]*zeta[1] *(b * b + d * d) <  dist2(b, d, t, v, d1, d2, d3))
+//		if ((max(b, t) <= zeta[2] * min(d, v))
+//				&& (min(dist_s2(b, t, d1), dist_s2(d, v, d2)) > 0))
 //			mexPrintf("Qx1x2");
-			return calcParIntQY1X1(b, d, t, v, d1, d2, d3);
+		return calcParIntQX1X2(b, d, t, v, d1, d2, d3);
 
-	if ((t * t + v * v) <= zeta[0] * dist2(b, d, t, v, d1, d2, d3))
+	if (zeta[0]*zeta[0] *(t * t + v * v) <  dist2(b, d, t, v, d1, d2, d3))
 		return calcParIntQ(b, d, t, v, d1, d2, d3);
 
 	return calcParIntA(b, d, t, v, d1, d2, d3);
diff --git a/src/test_sol.m b/src/test_sol.m
index e1d29e2..598537b 100644
--- a/src/test_sol.m
+++ b/src/test_sol.m
@@ -9,11 +9,11 @@ mex mex_build_V.cpp slpRectangle.cpp CFLAGS="\$CFLAGS -fopenmp" CXXFLAGS="\$CXXF
 % Test ausführen
 
 %Anzahl der Schritte
-steps = 4;
+steps = 32;
 
 %Art der Berechnungen
-type = [2 2 2 2 2];
-zeta = {[0 2] [1 2] [2 2] [3 2] [4 2]};
+type = [1 3];
+zeta = { [2 2 2] [2 2 2]};
 
 
 %% Quad adaptiv anisotrop