[src] begonnen mit Suche nach eigenartiger Verfeinerung

[src] ein paar Kommentare bearbeitet
[src] angefangen refineQuad noch mal zu überarbeiten -> ZIEL: nur verfeinern wenn wirklich nötig

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

diff --git a/src/A_step.m b/src/A_step.m
index 0b20dd44a94fe8efd910125518909662a63f6114..a76b4c3dc08a6c28a9b9c792bfb4dc19123f9f46 100644 (file)
--- a/src/A_step.m
+++ b/src/A_step.m
@@ -25,7 +25,7 @@ time = zeros(1,3);
   [coordinates_fine,elements_fine,neigh_fine,f2s]=refineQuad(G_C,G_E,G_N,2);
   time(1) = toc;
 
-  b = sqrt(sum(quadNorm(coordinates_fine,elements_fine,'w').^2,2));
+  b = sqrt(sum(quadNorm(coordinates_fine,elements_fine,'w').^2,2))
 
   tic
   data = size(G_E,1);
@@ -50,8 +50,15 @@ time = zeros(1,3);
   end
   time(2) = toc;
   
+  % nur RandElemente Verfeinern
+%   marked = ones(1,size(G_E,1));
+%   marked(find(sum((G_N(:,1:4)==0),2))) = 2;
   
-  marked = mark(x_fine(f2s)',ind,eta,eps);
+  
+  % Adaptiv Verfeinern
+  marked = mark(x_fine(f2s)',ind,eta,eps)
+%   plotShape(coordinates_fine,elements_fine,x_fine);
+%  plotMark(find(marked>1),G_C,G_E);
   
 %   clear 'coordinates_fine' 'elements_fine' 'neigh_fine' 'f2s'
   

diff --git a/src/computeEstSlpMuTilde.m b/src/computeEstSlpMuTilde.m
index 76c6f5738ff49ae37f81ed72f42111e5a69c7c59..13af47a32ade65cfca43c8d44437b9f8a9efe884 100644 (file)
--- a/src/computeEstSlpMuTilde.m
+++ b/src/computeEstSlpMuTilde.m
@@ -1,4 +1,10 @@
 function ind=computeEstSlpMuTilde(x_fine,coo,ele,f2s)
+% x_fine - Phi
+% coo - coordinaten
+% ele - elemente
+% f2s - Vater Sohn beziehungen
+
+
 xF2S = x_fine(f2s)';
 
 if(size(xF2S,1)==1)

diff --git a/src/mark.m b/src/mark.m
index 83af46fa1927d73c8bbc781fac3a3f480ed6e848..01d240eadfb22ca4c7521dac52a0c0518f868962 100644 (file)
--- a/src/mark.m
+++ b/src/mark.m
@@ -1,9 +1,10 @@
 function REF = mark(xF2S,ind,theta,eta)
 % function REF = mark(xF2S,ind,theta,eta)
-% xF2S - Father son 
+% xF2S - Father son relation
 % ind - error estimator
 % eta - refine element? (0..1, 1 = All)
 % eps - refine how? (0 = Isotrop)
+% REF - vector with entries [1 : 4] 
 
 if(size(xF2S,1)==1)
   xF2S = xF2S';
@@ -25,18 +26,18 @@ if(theta <1)
 
     ell = find(sum_ind  >= sum_ind(end) * theta,1);
 
-    t1(idx(ell:end)) = 1;
+    t1(idx(ell:end)) = 1;   % Nicht verfeinert
 end
 
 
 %% Wie muss verfeinert werden
-if(eps > 0)
+if(eps > 0) % Horizontal oder Vertikal
     t3 = (eta*abs(Ct(3,:)) >= sqrt(Ct(2,:).^2+Ct(4,:).^2));
     REF(t3) = 3;
     t4 = (eta*abs(Ct(4,:)) >= sqrt(Ct(2,:).^2+Ct(3,:).^2));
     REF(t4) = 4;
 end
-REF(~(t4+t3+t1)) = 2;
+REF(~(t4+t3+t1)) = 2;   % Rest wird Horizontal UND Vertikal geteilt
 
 
 end

diff --git a/src/plotShape.m b/src/plotShape.m
index cf5f42e699eb10b7cc0b1ae9476d07606f99ce5e..d50354d5d003fb63901acb0630f55b4c64dd7769 100644 (file)
--- a/src/plotShape.m
+++ b/src/plotShape.m
@@ -72,6 +72,15 @@ if(n)
     end
 end
 
+%     anorm = quadNorm(coordinates,elements);
+%     for idx = 1:eles
+%         current = sum(coordinates(elements(idx,[2,4])',:),1)/2;
+%         current = [current ; current+anorm(idx,:)*ind(idx);coordinates(elements(idx,1)',:)];
+%         plot3(current(:,1),current(:,2),current(:,3),'r'); % Zeichnet Oberflaeche
+%         scatter3(current(2,1),current(2,2),current(2,3),'xr');
+%         hold on
+%     end
+
 xlabel 'x'
 ylabel 'y'
 zlabel 'z'

diff --git a/src/refineQuad.m b/src/refineQuad.m
index 9a2f7d9c11f358b1d1f02c1f23f111d966f24005..447051448a500d394e8fd923ed836a7c1d182323 100644 (file)
--- a/src/refineQuad.m
+++ b/src/refineQuad.m
@@ -71,12 +71,13 @@ while(1==1)
     break;
   end
   for ele = ref(randperm(length(ref)))  %LOL Zufall bringts
+      
  % # HangingNode Check
     N = G_ref_N(ele,G_ref_N(ele,5:8)==0);
     N2 = N(N~=0); %Nachbarn der Kanten mit nur einem Nachbar
     if(~isempty(N2))
       N3=mod(find((G_ref_N(N2',:)==ele)')-1,4)+5; % ACHTUNG noch mal überprüfen
-      N4=N2(diag(G_ref_N(N2',N3'))~=0);
+      N4=N2(diag(G_ref_N(N2',N3'))~=0)
         if(~isempty(N4))
           for i = N4
             if(i>length(G_ref_t))
@@ -255,10 +256,6 @@ global G_ref_E;
 global G_ref_N;
 global G_ref_f2s;
 
-
-%Innere Beziehungen setzen
-% ...
-
 % plotShape(G_ref_C,G_ref_E);
 
 this = G_ref_N(ele,:);

diff --git a/src/slpRectangle.cpp b/src/slpRectangle.cpp
index a6d212cc7ef5d64b9e05cfb9c860aaa4b973ec8a..79a2dde10cf547e8e2426e71f5d6b8e78d46db10 100644 (file)
--- a/src/slpRectangle.cpp
+++ b/src/slpRectangle.cpp
@@ -9,12 +9,9 @@
 #include "slpRectangle.hpp"\r
 #include "gauss.hpp"\r
 \r
-#define quad 2 // Anzahl der Quadratur Auswertungstellen 2^quad\r
-using namespace std;\r
-//using namespace slpR;\r
+#define quad 4 // Anzahl der Quadratur Auswertungstellen 2^quad\r
 \r
-//double gauss_w5[] = { 0.1185, 0.2393, 0.2844, 0.2393, 0.1185 };\r
-//double gauss_n5[] = { 0.0469, 0.2308, 0.5000, 0.7692, 0.9531 };\r
+using namespace std;\r
 \r
 int sign(double);\r
 //double arsinh(double);\r
@@ -42,7 +39,7 @@ double inline f_A(double x1, double x2, double y1, double y2, double l) {
        return 1. / sqrt((x1 - y1) * (x1 - y1) + (x2 - y2) * (x2 - y2) + l * l);\r
 }\r
 \r
-double slpR::g_QY(double p, double y, double x, double l) {\r
+double inline g_QY(double p, double y, double x, double l) {\r
 \r
        double sol = 0;\r
 \r
@@ -57,7 +54,7 @@ double slpR::g_QY(double p, double y, double x, double l) {
 }\r
 \r
 //y-x muss != 0 sein\r
-double slpR::g_AY(double p, double y, double x, double l) {\r
+double inline g_AY(double p, double y, double x, double l) {\r
        //printf("%.1f | %.1f | %.1f | %.1f +",p,x,y,l);\r
 \r
        double sol = 0;\r
@@ -88,7 +85,7 @@ double slpR::g_AY(double p, double y, double x, double l) {
        return sol;\r
 }\r
 \r
-double slpR::G_QY1Y2(double p, double y1, double y2, double x1, double x2,\r
+double inline G_QY1Y2(double p, double y1, double y2, double x1, double x2,\r
                double l) {\r
 \r
        double sol = 0;\r
@@ -106,7 +103,7 @@ double slpR::G_QY1Y2(double p, double y1, double y2, double x1, double x2,
        return sol;\r
 }\r
 \r
-double slpR::G_AY2X2(double y1, double y2, double x1, double x2, double d3) {\r
+double inline G_AY2X2(double y1, double y2, double x1, double x2, double d3) {\r
 \r
 //     cout << "("     << y1 << "," << y2 << ")(" << x1 << "," << x2 << ")(" << d3 << ")" << endl;\r
 \r
@@ -127,7 +124,7 @@ double slpR::G_AY2X2(double y1, double y2, double x1, double x2, double d3) {
        return sol;\r
 }\r
 \r
-double slpR::G_AY1Y2(double p, double y1, double y2, double x1, double x2,\r
+double inline G_AY1Y2(double p, double y1, double y2, double x1, double x2,\r
                double l) {\r
        //      printf("%.1f | %.1f %.1f | %.1f %.1f | %.1f +",p,x1,x2,y1,y2,l);\r
        double pt = p + 1.5;\r
@@ -164,14 +161,14 @@ double slpR::G_AY1Y2(double p, double y1, double y2, double x1, double x2,
        return sol;\r
 }\r
 \r
-double slpR::Gs_AX2Y1Y2(double p, double y1, double y2, double x1, double x2,\r
+double inline Gs_AX2Y1Y2(double p, double y1, double y2, double x1, double x2,\r
                double l) {\r
        double sol = 0;\r
 \r
        return sol;\r
 }\r
 \r
-double slpR::F_par(double x1, double x2, double y1, double y2, double d1,\r
+double inline F_par(double x1, double x2, double y1, double y2, double d1,\r
                double d2, double d3) {\r
 \r
        //       printf("%.1f %.1f | %.1f %.1f | %.1f %.1f %.1f",x1,x2,y1,y2,d1,d2,d3);\r
@@ -196,7 +193,7 @@ double slpR::F_par(double x1, double x2, double y1, double y2, double d1,
        return sol;\r
 }\r
 \r
-double slpR::F_ort(double x1, double x2, double y2, double y3, double d1,\r
+double inline F_ort(double x1, double x2, double y2, double y3, double d1,\r
                double d2, double d3) {\r
 \r
        //       printf("%.1f %.1f | %.1f %.1f | %.1f %.1f %.1f",x1,x2,y2,y3,d1,d2,d3);\r
@@ -248,7 +245,7 @@ double slpR::F_ort(double x1, double x2, double y2, double y3, double d1,
 \r
  }*/\r
 \r
-double slpR::apply0Int4(\r
+double inline apply0Int4(\r
                double(*f)(double, double, double, double, double, double, double),\r
                double b, double d, double t, double v, double d1, double d2,\r
                double d3) {\r
@@ -264,7 +261,7 @@ double slpR::apply0Int4(
 \r
 }\r
 \r
-double slpR::apply0Int2(\r
+double inline apply0Int2(\r
                double(*f)(double, double, double, double, double, double, double),\r
                double b, double d, double t, double v, double d1, double d2,\r
                double d3) {\r

diff --git a/src/slpRectangle.hpp b/src/slpRectangle.hpp
index fae183d408b9f1b317846c36debc92147ea81741..b7ff93567a4b4f6d768d9f7b87163898b4a2afb9 100644 (file)
--- a/src/slpRectangle.hpp
+++ b/src/slpRectangle.hpp
@@ -4,30 +4,30 @@
 namespace slpR {
 
        // sol = g0(p,y,x,l);
-       double g_AY(double, double, double, double);
+//     double g_AY(double, double, double, double);
        // sol = g0(p,y,x,l);
-       double g_QY(double, double, double, double);
+//     double g_QY(double, double, double, double);
 
        // sol = G00(p,y1,y2,x1,x2,l);
-       double G_AY1Y2(double, double, double, double, double, double);
+//     double G_AY1Y2(double, double, double, double, double, double);
        // sol = G00(p,y1,y2,x1,x2,l);
-       double G_AY2X2(double, double, double, double, double);
+//     double G_AY2X2(double, double, double, double, double);
        // sol = G00(p,y1,y2,x1,x2,l);
-       double G_QY1Y2(double, double, double, double, double, double);
-       double Gs_AX2Y1Y2(double, double, double, double, double, double);
+//     double G_QY1Y2(double, double, double, double, double, double);
+//     double Gs_AX2Y1Y2(double, double, double, double, double, double);
 
        // sol = F_par(x1,x2,y1,y2,d1,d2,d3);
-       double F_par(double, double, double, double, double, double, double);
+//     double F_par(double, double, double, double, double, double, double);
        // sol = F_ort(x1,x2,y2,y3,d1,d2,d3);
-       double F_ort(double, double, double, double, double, double, double);
+//     double F_ort(double, double, double, double, double, double, double);
 
        // sol = quad0Int4((F_par/F_ort),b,d,t,v,d1,d2,d3);
-       double apply0Int4(
-                       double(*)(double, double, double, double, double, double, double),
-                       double, double, double, double, double, double, double);
-       double apply0Int2(
-                       double(*)(double, double, double, double, double, double, double),
-                       double, double, double, double, double, double, double);
+//     double apply0Int4(
+//double(*)(double, double, double, double, double, double, double),
+       //double, double, double, double, double, double, double);
+                       //      double apply0Int2(
+                       //double(*)(double, double, double, double, double, double, double),
+                       //double, double, double, double, double, double, double);
 
        // sol = calcParInt.(b,d,t,v,d1,d2,d3);
        double calcParIntA(double, double, double, double, double, double, double);
@@ -40,6 +40,8 @@ namespace slpR {
        double calcOrtIntA(double, double, double, double, double, double, double);
        double calcOrtIntQX1X2(double, double, double, double, double, double, double);
 
+
+
        //Voll Analytische Berechnung
        double calcParIntO0(double, double, double, double, double, double, double,
                        double);