From: treecity <treecity@26120e32-c555-405d-b3e1-1f783fb42516>
Date: Mon, 2 May 2011 08:47:40 +0000 (+0000)
Subject: [build_A] weitere Fehler behoben
X-Git-Url: https://git.leopard-lacewing.eu/?a=commitdiff_plain;h=7255c5fd57aa95ea0bc529438726416429769027;p=bacc.git

[build_A] weitere Fehler behoben
[src] solveError weiter gebaut... Funktioniert aber immer noch nicht

git-svn-id: https://drops.fb12.tu-berlin.de/svn/bacc/trunk@17 26120e32-c555-405d-b3e1-1f783fb42516
---

diff --git a/src/build_A.cpp b/src/build_A.cpp
index 5896472..ada1f31 100644
--- a/src/build_A.cpp
+++ b/src/build_A.cpp
@@ -132,7 +132,8 @@ for (int j=0;j<em;++j){
 	 if(rx=ry){
 		 if(rx=2){
 			 //printf("%d%d %.1f %.1f | %.1f %.1f || %.1f %.1f | %.1f %.1f || %.1f %.1f %.1f\n",j,k,x1[0],x1[1],x2[0],x3[1],y1[0],y1[1],y2[0],y3[1],d[0],d[1],d[2]);
-			 tmp = quadInt(F_par,x1[0],x1[1],x2[0],x3[1],y1[0],y1[1],y2[0],y3[1],d[0],d[1],d[2]);
+			 //printf("%d%d %.1f %.1f %.1f\n",j,k,d[0],d[1],d[2]);
+             tmp = quadInt(F_par,x1[0],x1[1],x2[0],x3[1],y1[0],y1[1],y2[0],y3[1],d[0],d[1],d[2]);
 		//	 printf("%d%d|%.2f\n",j,k,tmp);
 			 A[(k*em)+j] = 1./(4*my_PI)*tmp;
 
@@ -184,6 +185,7 @@ for (int j=0;j<em;++j){
         return log(x+sqrt(x*x+1));
     }
     
+    //y-x muss != 0 sein
     double g0(double p, double y, double x, double l){
     //printf("%.1f | %.1f | %.1f | %.1f +",p,x,y,l);
     
@@ -191,8 +193,8 @@ for (int j=0;j<em;++j){
 
 		if(l!=0){
             if(p==0.5){
-                sol = (y-x)/2*sqrt((y-x)*(y-x)+
-                        l*l)+l*l/2*arsinh((y-x)/abs(l));
+                sol = (y-x)/2*sqrt((y-x)*(y-x)+l*l)+
+                        l*l/2*arsinh((y-x)/abs(l));
       //          printf("%.2f |",sol);
             }else if(p==0)
                 sol = y-x;
@@ -203,8 +205,8 @@ for (int j=0;j<em;++j){
             else if(p==-1.5)
                 sol = (y-x)/((l*l)*sqrt((y-x)*(y-x)+l*l));
             else
-                sol = (y-x)*pow((y-x)*(y-x)+l*l,p)+2*p*l*
-                        g0(p-1,y,x,l);
+                sol = (y-x)*pow((y-x)*(y-x)+l*l,p)+2*p*l*l*
+                        g0(p-1,y,x,l)/(2*p+1);
 		}else{
             if(p==-0.5)
 				sol = sign(y-x)*log(abs(y-x));
@@ -261,7 +263,7 @@ for (int j=0;j<em;++j){
             sol -= (x2-y2-d2)*g0(0.5,x2,y2+d2,sqrt((x1-y1-d1)*(x1-y1-d1)+d3*d3));
 
         double hlp = ((x1-y1-d1)*(x1-y1-d1)+(x2-y2-d2)*(x2-y2-d2)+d3*d3);
-        sol +=  1/3.*hlp*sqrt(hlp);
+        sol +=  1./3*hlp*sqrt(hlp);
         return sol;
     }
     
diff --git a/src/test_solve1.mat b/src/test_solve1.mat
deleted file mode 100644
index d220296..0000000
Binary files a/src/test_solve1.mat and /dev/null differ
diff --git a/src/test_solveError.m b/src/test_solveError.m
index b3ee998..848a706 100644
--- a/src/test_solveError.m
+++ b/src/test_solveError.m
@@ -1,25 +1,44 @@
 
-load exmpl_2DLShape
-
-re_times = 6;
-r_type = 1;
-
-xNorm = zeros(re_times+1,2);
-
-    test_solve
-    xNorm(1,:) = [size(x,1) xe];
 
-for l = 2:re_times
-    test_refine
-    
-    test_solve
-    xNorm(l,:) = [size(x,1) xe];
-end
-
-
-loglog(xNorm(:,1),xNorm(:,2));
+load exmpl_2DLShape
 
-legend('x^TAx')
-xlabel('Elemente')
-ylabel('xNorm')
+A = build_A(coordinates,elements);
+b = sqrt(sum(quad_norm(coordinates,elements,'w').^2,2));
+x = A\b;
+xe = x'*A*x;
+
+[coordinates_fine,elements_fine, f2s]=refineQuad(coordinates,elements,ones(1,size(elements,1)));
+A_fine = build_A(coordinates_fine,elements_fine);
+b = sqrt(sum(quad_norm(coordinates_fine,elements_fine,'w').^2,2));
+x_fine = A_fine\b;
+xe_fine = x_fine'*A_fine*x_fine;
+
+l = x_fine(f2s(1,:));
+T4 = [1 1 1 1; 1 -1 1 -1; 1 1 -1 -1;1 -1 -1 1];
+l
+T4*l*1/4
+
+[xe xe_fine]
+
+% re_times = 6;
+% r_type = 1;
+% 
+% xNorm = zeros(re_times+1,2);
+% 
+%     test_solve
+%     xNorm(1,:) = [size(x,1) xe];
+% 
+% for l = 2:re_times
+%     
+%     
+%     test_solve
+%     xNorm(l,:) = [size(x,1) xe];
+% end
+% 
+% 
+% loglog(xNorm(:,1),xNorm(:,2));
+% 
+% legend('x^TAx')
+% xlabel('Elemente')
+% ylabel('xNorm')