[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 58964728d19453af22925f3ce23d9573fe7c9d74..ada1f31336c74662b2c190c42ee2d1d7448d4fa0 100644 (file)
--- a/src/build_A.cpp
+++ b/src/build_A.cpp
@@ -132,7 +132,8 @@ for (int j=0;j<em;++j){
         if(rx=ry){\r
                 if(rx=2){\r
                         //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]);\r
-                        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]);\r
+                        //printf("%d%d %.1f %.1f %.1f\n",j,k,d[0],d[1],d[2]);\r
+             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]);\r
                //       printf("%d%d|%.2f\n",j,k,tmp);\r
                         A[(k*em)+j] = 1./(4*my_PI)*tmp;\r
 \r
@@ -184,6 +185,7 @@ for (int j=0;j<em;++j){
         return log(x+sqrt(x*x+1));\r
     }\r
     \r
+    //y-x muss != 0 sein\r
     double g0(double p, double y, double x, double l){\r
     //printf("%.1f | %.1f | %.1f | %.1f +",p,x,y,l);\r
     \r
@@ -191,8 +193,8 @@ for (int j=0;j<em;++j){
 \r
                if(l!=0){\r
             if(p==0.5){\r
-                sol = (y-x)/2*sqrt((y-x)*(y-x)+\r
-                        l*l)+l*l/2*arsinh((y-x)/abs(l));\r
+                sol = (y-x)/2*sqrt((y-x)*(y-x)+l*l)+\r
+                        l*l/2*arsinh((y-x)/abs(l));\r
       //          printf("%.2f |",sol);\r
             }else if(p==0)\r
                 sol = y-x;\r
@@ -203,8 +205,8 @@ for (int j=0;j<em;++j){
             else if(p==-1.5)\r
                 sol = (y-x)/((l*l)*sqrt((y-x)*(y-x)+l*l));\r
             else\r
-                sol = (y-x)*pow((y-x)*(y-x)+l*l,p)+2*p*l*\r
-                        g0(p-1,y,x,l);\r
+                sol = (y-x)*pow((y-x)*(y-x)+l*l,p)+2*p*l*l*\r
+                        g0(p-1,y,x,l)/(2*p+1);\r
                }else{\r
             if(p==-0.5)\r
                                sol = sign(y-x)*log(abs(y-x));\r
@@ -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));\r
 \r
         double hlp = ((x1-y1-d1)*(x1-y1-d1)+(x2-y2-d2)*(x2-y2-d2)+d3*d3);\r
-        sol +=  1/3.*hlp*sqrt(hlp);\r
+        sol +=  1./3*hlp*sqrt(hlp);\r
         return sol;\r
     }\r
     \r

diff --git a/src/test_solve1.mat b/src/test_solve1.mat
deleted file mode 100644 (file)
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 b3ee99865ec128a2767729a96777730c3845f903..848a7069b54b63e8c0042a16686be0d6bb108eed 100644 (file)
--- a/src/test_solveError.m
+++ b/src/test_solveError.m
@@ -1,25 +1,44 @@
 \r
-load exmpl_2DLShape\r
-\r
-re_times = 6;\r
-r_type = 1;\r
-\r
-xNorm = zeros(re_times+1,2);\r
-\r
-    test_solve\r
-    xNorm(1,:) = [size(x,1) xe];\r
 \r
-for l = 2:re_times\r
-    test_refine\r
-    \r
-    test_solve\r
-    xNorm(l,:) = [size(x,1) xe];\r
-end\r
-\r
-\r
-loglog(xNorm(:,1),xNorm(:,2));\r
+load exmpl_2DLShape\r
 \r
-legend('x^TAx')\r
-xlabel('Elemente')\r
-ylabel('xNorm')\r
+A = build_A(coordinates,elements);\r
+b = sqrt(sum(quad_norm(coordinates,elements,'w').^2,2));\r
+x = A\b;\r
+xe = x'*A*x;\r
+\r
+[coordinates_fine,elements_fine, f2s]=refineQuad(coordinates,elements,ones(1,size(elements,1)));\r
+A_fine = build_A(coordinates_fine,elements_fine);\r
+b = sqrt(sum(quad_norm(coordinates_fine,elements_fine,'w').^2,2));\r
+x_fine = A_fine\b;\r
+xe_fine = x_fine'*A_fine*x_fine;\r
+\r
+l = x_fine(f2s(1,:));\r
+T4 = [1 1 1 1; 1 -1 1 -1; 1 1 -1 -1;1 -1 -1 1];\r
+l\r
+T4*l*1/4\r
+\r
+[xe xe_fine]\r
+\r
+% re_times = 6;\r
+% r_type = 1;\r
+% \r
+% xNorm = zeros(re_times+1,2);\r
+% \r
+%     test_solve\r
+%     xNorm(1,:) = [size(x,1) xe];\r
+% \r
+% for l = 2:re_times\r
+%     \r
+%     \r
+%     test_solve\r
+%     xNorm(l,:) = [size(x,1) xe];\r
+% end\r
+% \r
+% \r
+% loglog(xNorm(:,1),xNorm(:,2));\r
+% \r
+% legend('x^TAx')\r
+% xlabel('Elemente')\r
+% ylabel('xNorm')\r
 \r