[src] eigenartig -> Funktionen anders füttern geht aber Variablen tauschen nicht...

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

diff --git a/src/mex_build_As1.cpp b/src/mex_build_As1.cpp
index d54e5126b103f17b7fe30ced6058419502bd5c26..485b74fb199f9a2ac2b4ac337c6412896cd8b5c0 100644 (file)
--- a/src/mex_build_As1.cpp
+++ b/src/mex_build_As1.cpp
@@ -69,7 +69,7 @@ void mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[]) {
        double * d = mxGetPr(mxCreateDoubleMatrix(3, 1, mxREAL));\r
        double * dt = mxGetPr(mxCreateDoubleMatrix(3, 1, mxREAL));\r
 \r
-       double tmp;\r
+       double tmp,tmp2;\r
 \r
        int itmp;\r
 \r
@@ -206,11 +206,11 @@ void mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[]) {
 \r
                                if(xa[rxa]*xa[rxa]+xb[rxb]*xb[rxb]>ya[rya]*ya[rya]+yb[ryb]*yb[ryb]){\r
 \r
-       /*                              printf("(%d,%d@%d,%d)",rxa,rxb,rya,ryb);\r
+                                       printf("+(%d,%d@%d,%d)",rxa,rxb,rya,ryb);\r
                                        printf("[%.2f,%.2f@%.2f,%.2f]",xa[rxa],xb[rxb],ya[rya],yb[ryb]);\r
                                        printf("{%d%d} %.1f %.1f %.1f",j,k,d[0],d[1],d[2]);\r
-                                       printf(" | ");\r
-*/\r
+                                       printf("\n");\r
+\r
 \r
                                        vtmp = xa; xa = ya; ya = vtmp;\r
                                        vtmp = xb; xb = yb; yb = vtmp;\r
@@ -222,12 +222,12 @@ void mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[]) {
                                        for (i = 0; i<3;++i)\r
                                                d[i] = -d[i];\r
 \r
-                                       /*\r
 \r
-                                       printf("(%d,%d@%d,%d)",rxa,rxb,rya,ryb);\r
+\r
+                                       printf("-(%d,%d@%d,%d)",rxa,rxb,rya,ryb);\r
                                        printf("[%.2f,%.2f@%.2f,%.2f]",xa[rxa],xb[rxb],ya[rya],yb[ryb]);\r
                                        printf("{%d%d} %.1f %.1f %.1f",j,k,d[0],d[1],d[2]);\r
-                                       printf("\n");*/\r
+                                       printf("\n");\r
                                }\r
 \r
 \r
@@ -243,10 +243,23 @@ void mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[]) {
                                        //printf("%d%d %.1f %.1f | %.1f %.1f || %.1f %.1f | %.1f %.1f || %.1f %.1f %.1f\n",j,k,x0[0], x1[0], x0[1], x3[1], y0[0],\r
                                        //              y1[0], y0[1], 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
+\r
+                                       if(xa[rxa]*xa[rxa]+xb[rxb]*xb[rxb]<=ya[rya]*ya[rya]+yb[ryb]*yb[ryb]){\r
+\r
                                        tmp\r
                                                        = apply0Int(F_par, fabs(xa[rxa]), fabs(xb[rxb]),\r
-                                                                       fabs(ya[rxa]), fabs(yb[rxb]), d[rxa],\r
+                                                                       fabs(ya[rya]), fabs(yb[ryb]), d[rxa],\r
                                                                        d[rxb], d[rx]);\r
+\r
+                                       }else{\r
+\r
+                                       tmp\r
+                                                       = apply0Int(F_par, fabs(ya[rya]), fabs(yb[ryb]),\r
+                                                                       fabs(xa[rxa]), fabs(xb[rxb]), -d[rxa],\r
+                                                                       -d[rxb], -d[rx]);\r
+                                       }\r
+                       //              if(fabs(tmp-tmp2)>10e-16)\r
+                       //                      printf("wtf");\r
                                        //       printf("%d%d|%.2f\n",j,k,tmp);\r
                                } else {\r
                                        //printf("%d%d %.1f %.1f | %.1f %.1f || %.1f %.1f | %.1f %.1f || %.1f %.1f %.1f\n",j,k,x0[0], x1[0], x0[1], x3[1], y0[0],\r

diff --git a/src/save1.mat b/src/save1.mat
new file mode 100644 (file)
index 0000000..89b1459
Binary files /dev/null and b/src/save1.mat differ

diff --git a/src/save10.mat b/src/save10.mat
new file mode 100644 (file)
index 0000000..d246641
Binary files /dev/null and b/src/save10.mat differ

diff --git a/src/save11.mat b/src/save11.mat
new file mode 100644 (file)
index 0000000..093cc3d
Binary files /dev/null and b/src/save11.mat differ

diff --git a/src/save12.mat b/src/save12.mat
new file mode 100644 (file)
index 0000000..095aa94
Binary files /dev/null and b/src/save12.mat differ

diff --git a/src/save13.mat b/src/save13.mat
new file mode 100644 (file)
index 0000000..e028544
Binary files /dev/null and b/src/save13.mat differ

diff --git a/src/save2.mat b/src/save2.mat
new file mode 100644 (file)
index 0000000..42820ed
Binary files /dev/null and b/src/save2.mat differ

diff --git a/src/save3.mat b/src/save3.mat
new file mode 100644 (file)
index 0000000..cade101
Binary files /dev/null and b/src/save3.mat differ

diff --git a/src/save4.mat b/src/save4.mat
new file mode 100644 (file)
index 0000000..2748262
Binary files /dev/null and b/src/save4.mat differ

diff --git a/src/save5.mat b/src/save5.mat
new file mode 100644 (file)
index 0000000..e4d3227
Binary files /dev/null and b/src/save5.mat differ

diff --git a/src/save6.mat b/src/save6.mat
new file mode 100644 (file)
index 0000000..1a71948
Binary files /dev/null and b/src/save6.mat differ

diff --git a/src/save7.mat b/src/save7.mat
new file mode 100644 (file)
index 0000000..8203440
Binary files /dev/null and b/src/save7.mat differ

diff --git a/src/save8.mat b/src/save8.mat
new file mode 100644 (file)
index 0000000..297ebba
Binary files /dev/null and b/src/save8.mat differ

diff --git a/src/save9.mat b/src/save9.mat
new file mode 100644 (file)
index 0000000..4ce9f2c
Binary files /dev/null and b/src/save9.mat differ

diff --git a/src/test_solveError_Grid.m b/src/test_solveError_Grid.m
new file mode 100755 (executable)
index 0000000..4c745b6
--- /dev/null
+++ b/src/test_solveError_Grid.m
@@ -0,0 +1,57 @@
+\r
+dataIso =[];\r
+dataAniso =[];\r
+maxSize = 500;\r
+Netz = 'exmpl_2DLShape'; %Energienorm sollte gegen 8.28466 konvergieren\r
+% Netz = 'exmpl_3DFichCube'; %Energienorm sollte gegen 162265 konvergieren\r
+% Netz = 'exmpl_2DQuad';\r
+\r
+%% Anisotrope Verfeinerung\r
+disp 'Anisotrop'\r
+load(Netz)\r
+\r
+ref = 0;\r
+while size(elements,1)<maxSize\r
+  tic\r
+  ref = ref+1; \r
+\r
+  %Netz Isotrop verfeinern (4 Teile)\r
+  [coordinates_fine,elements_fine, f2s]=refineQuad(coordinates,elements,2*ones(1,size(elements,1)));\r
+  \r
+  %Matrix mit feinem Netz berechnen\r
+  A_fine = mex_build_A(coordinates_fine,elements_fine);\r
+  \r
+  %rechte Seite Aufbauen\r
+  b = sqrt(sum(quadNorm(coordinates_fine,elements_fine,'w').^2,2));\r
+  \r
+  %Lösungsvektor bestimmen\r
+  x_fine = A_fine\b;\r
+  \r
+  %Energienorm^2 bestimmen\r
+  xe_fine = x_fine'*A_fine*x_fine;\r
+\r
+  % Welche Elemente sollen verfeinert werden\r
+  ind = computeEstSlpMuTilde(x_fine,coordinates,elements,f2s);\r
+  \r
+  % Wie soll verfeinert werden\r
+  marked = mark(x_fine(f2s)',ind,0.5,0.5);\r
+\r
+  %Daten zur Auswertung zwischenspeichern (testweise?)\r
+  dataAniso(ref,:) = [size(elements,1) sqrt(sum(ind)) xe_fine];\r
+  \r
+  %Sicherheit falls keine Elemente markiert sind\r
+  if(sum(marked~=1)==0)\r
+    disp 'Abgebrochen weil nicht mehr verfeinert wird!'\r
+    break;\r
+  end\r
+  \r
+  %Verfeinern\r
+  [coordinates, elements] = refineQuad(coordinates,elements,marked);\r
+  file = ['save' int2str(ref)];\r
+  save ([file], 'coordinates', 'elements')\r
+  toc\r
+end\r
+dataAniso\r
+\r
+clear ref b A_fine coordinates_fine elements_fine marked ind f2s x_fine\r
+\r

diff --git a/src/test_solveError_Grid2.m b/src/test_solveError_Grid2.m
new file mode 100755 (executable)
index 0000000..f90d5af
--- /dev/null
+++ b/src/test_solveError_Grid2.m
@@ -0,0 +1,126 @@
+\r
+dataIso =[];\r
+dataAniso =[];\r
+maxSize = 500;\r
+fileName = './test_save';\r
+Netz = 'exmpl_2DLShape'; %Energienorm sollte gegen 8.28466 konvergieren\r
+% Netz = 'exmpl_3DFichCube'; %Energienorm sollte gegen 162265 konvergieren\r
+% Netz = 'exmpl_2DQuad';\r
+\r
+\r
+% %% Anisotrope Verfeinerung\r
+% disp 'Anisotrop1'\r
+% load(Netz)\r
+% \r
+% ref = 0;\r
+% while size(elements,1)<maxSize\r
+%   tic\r
+%   ref = ref+1; \r
+%   \r
+% \r
+% \r
+%   %Netz Isotrop verfeinern (4 Teile)\r
+%   [coordinates_fine,elements_fine, f2s]=refineQuad(coordinates,elements,2*ones(1,size(elements,1)));\r
+%   \r
+%   %Matrix mit feinem Netz berechnen\r
+%   A_fine = mex_build_A(coordinates_fine,elements_fine);\r
+%   \r
+%   %rechte Seite Aufbauen\r
+%   b = sqrt(sum(quadNorm(coordinates_fine,elements_fine,'w').^2,2));\r
+%   \r
+%   %Lösungsvektor bestimmen\r
+%   x_fine = A_fine\b;\r
+%   \r
+%   %Energienorm^2 bestimmen\r
+%   xe_fine = x_fine'*A_fine*x_fine;\r
+% \r
+%   % Welche Elemente sollen verfeinert werden\r
+%   ind = computeEstSlpMuTilde(x_fine,coordinates,elements,f2s);\r
+%   \r
+%   % Wie soll verfeinert werden\r
+%   marked = mark(x_fine(f2s)',ind,0.5,0.5);\r
+% \r
+%   %Daten zur Auswertung zwischenspeichern (testweise?)\r
+%   dataAniso(ref,:) = [size(elements,1) sqrt(sum(ind)) xe_fine];\r
+%   \r
+%   %Sicherheit falls keine Elemente markiert sind\r
+%   if(sum(marked~=1)==0)\r
+%     disp 'Abgebrochen weil nicht mehr verfeinert wird!'\r
+%     break;\r
+%   end\r
+%   \r
+%   %Verfeinern\r
+%   %[coordinates, elements] = refineQuad(coordinates,elements,marked);\r
+%   file = ['save' int2str(ref)];\r
+%   load(file);\r
+%   toc\r
+% end\r
+% d1 = dataAniso\r
+\r
+\r
+%% Anisotrope Verfeinerung\r
+disp 'Anisotrop2'\r
+load(Netz)\r
+\r
+ref = 0;\r
+while size(elements,1)<maxSize\r
+  tic\r
+  ref = ref+1; \r
+  \r
+\r
+\r
+  %Netz Isotrop verfeinern (4 Teile)\r
+  [coordinates_fine,elements_fine, f2s]=refineQuad(coordinates,elements,2*ones(1,size(elements,1)));\r
+  \r
+  %Matrix mit feinem Netz berechnen\r
+  A_fine = mex_build_As1(coordinates_fine,elements_fine);\r
+  \r
+  %rechte Seite Aufbauen\r
+  b = sqrt(sum(quadNorm(coordinates_fine,elements_fine,'w').^2,2));\r
+  \r
+  %Lösungsvektor bestimmen\r
+  x_fine = A_fine\b;\r
+  \r
+  %Energienorm^2 bestimmen\r
+  xe_fine = x_fine'*A_fine*x_fine;\r
+\r
+  % Welche Elemente sollen verfeinert werden\r
+  ind = computeEstSlpMuTilde(x_fine,coordinates,elements,f2s);\r
+  \r
+  % Wie soll verfeinert werden\r
+  marked = mark(x_fine(f2s)',ind,0.5,0.5);\r
+\r
+  %Daten zur Auswertung zwischenspeichern (testweise?)\r
+  dataAniso(ref,:) = [size(elements,1) sqrt(sum(ind)) xe_fine];\r
+  \r
+  %Sicherheit falls keine Elemente markiert sind\r
+  if(sum(marked~=1)==0)\r
+    disp 'Abgebrochen weil nicht mehr verfeinert wird!'\r
+    break;\r
+  end\r
+  \r
+  %Verfeinern\r
+  %[coordinates, elements] = refineQuad(coordinates,elements,marked);\r
+  file = ['save' int2str(ref)];\r
+  load(file);\r
+  toc\r
+end\r
+d2 = dataAniso\r
+\r
+clear ref b A_fine coordinates_fine elements_fine marked ind f2s x_fine\r
+\r
+\r
+%% Alles Zeichnen\r
+figure (1)\r
+loglog(d1(:,1),d1(:,2),d2(:,1),d2(:,2));\r
+\r
+legend('Isotrop','Anisotrop')\r
+xlabel('Elemente')\r
+ylabel('Fehler')\r
+\r
+figure (2)\r
+loglog(d1(:,1),d1(:,3),d2(:,1),d2(:,3));\r
+\r
+legend('Isotrop','Anisotrop')\r
+xlabel('Elemente')\r
+ylabel('EnergieNorm')\r