From acc0373663b5090c3cc2d2f2bcb40ed59cab588b Mon Sep 17 00:00:00 2001 From: treecity Date: Sun, 10 Apr 2011 22:01:25 +0000 Subject: [PATCH] [src] g G und beide F in Matlab umgesetzt [src] kleine korrektur an der refineQuad (TypReihenfolge) git-svn-id: https://drops.fb12.tu-berlin.de/svn/bacc/trunk@6 26120e32-c555-405d-b3e1-1f783fb42516 --- src/F_ort.m | 14 ++++++++++++++ src/F_par.m | 11 +++++++++++ src/G00.m | 18 ++++++++++++++++++ src/g0.m | 31 +++++++++++++++++++++++++++++++ src/plotShape.m | 21 ++++++++++++++------- src/quad_norm.m | 5 +++-- src/refineQuad.m | 36 ++++++++++++++++++------------------ src/test_refine.m | 2 +- 8 files changed, 110 insertions(+), 28 deletions(-) create mode 100644 src/F_ort.m create mode 100644 src/F_par.m create mode 100644 src/G00.m create mode 100644 src/g0.m diff --git a/src/F_ort.m b/src/F_ort.m new file mode 100644 index 0000000..2c5636b --- /dev/null +++ b/src/F_ort.m @@ -0,0 +1,14 @@ +function sol = F_ort(x, y, delta) + + xd = x(1) -delta(1); + xyd=x(2)-y(1)-delta(2); + yd = y(2)+delta(3); + + sol = -G00(0.5,[y(2) x(1)]',[-delta(3) delta(1)],xyd)... + -xd*xyd*G00(-0.5,[x(2) y(2)],[yd -delta(3)],xd)... + +xd*g0(0.5,y(2),-delta(3),norm([xd xyd]))... + -yd*xyd*G00(-0.5,[x(1) x(2)],[delta(1) yd],-yd)... + +yd*g0(0.5,x(1),delta(1),norm([xyd yd])); + +end + diff --git a/src/F_par.m b/src/F_par.m new file mode 100644 index 0000000..c214e4d --- /dev/null +++ b/src/F_par.m @@ -0,0 +1,11 @@ +function sol = F_par(x, y, delta) + + xyd=x-y-delta(1:2); + + sol = xyd(1)*xyd(2)*G(-0.5,x,y+delta(1:2),delta(3)... + -xyd(1)*g(0.5,x(1),y(1)+delta(1),norm([xyd(2) delta(3)]))... + -xyd(2)*g(0.5,x(2),y(2)+delta(2),norm([xyd(1) delta(3)]))... + +1/3*(xyd(1)^2+xyd(2)^2+delta(3)^2)^(1.5)); + +end + diff --git a/src/G00.m b/src/G00.m new file mode 100644 index 0000000..777e6f5 --- /dev/null +++ b/src/G00.m @@ -0,0 +1,18 @@ +function sol = G00(p,y,x,lambda) + + sol = nan; + yx = y-x; + if(p==-1.5) + if(lambda==0) + sol = -norm(yx)/prod(yx); + else + sol = sign(prod(yx))/(2*abs(lambda))*... + acos(-2*prod(yx)^2/prod((yx).^2+lambda^2)+1); + end + elseif(p==-0.5) + sol = (2*p*lambda^2*G(p-1,y,x,lambda)... + + sum( yx(1).*g(p,y(2),x(2),norm([xy(1) lambda])))... + + sum( yx(2).*g(p,y(1),x(1),norm([xy(2) lambda]))))/(2*p+2); + end + +end diff --git a/src/g0.m b/src/g0.m new file mode 100644 index 0000000..c8e1632 --- /dev/null +++ b/src/g0.m @@ -0,0 +1,31 @@ +function sol = g0(p, y, x, lambda) + + sol = nan; + yx = y-x; + if(lambda~=0) + if(p==0.5) + sol = yx/2*((yx^2 + lambda^2)^0.5)+lambda^2/2*arsinh(yx/abs(lambda)); + elseif(p==0) + sol = yx; + elseif(p==-0.5) + sol = arsinh(yx/abs(lambda)); + elseif(p==-1) + sol = atan(yx/abs(lambda)); + elseif(p==-1.5) + sol = yx/(lambda^2)*((yx^2 + lambda^2)^-0.5); + else + sol = yx*((yx^2 + lambda^2)^p)+2*p*lambda^2*g(p-1,y,x,lambda); + end + else + if(p==-0.5) + sol = sign(yx)*log(abs(yx)); + else + sol = (yx).*abs(yx).^(2*p)/(2*p+1); + end + end +end + +function sol = arsinh(x) %inline + sol = log(x+norm([x 1])); +end + diff --git a/src/plotShape.m b/src/plotShape.m index d7c9a9a..8aedde4 100644 --- a/src/plotShape.m +++ b/src/plotShape.m @@ -3,11 +3,12 @@ function plotShape(coordinates, elements, varargin) % plotShape(coordinates,elements) % plotShape(coordinates,elements,'FLAG') % -% Diese Funktion Zeichnet alle Vierecke mit ausgef llten F l c h e n , wenn als -% dritter Parameter noch 'f' angegeben wird. Werden nur die Kanten -% gezeichnet. Wenn 'c' angegeben wird werden auch die Ecken als kleine rote -% Kreise dargestellt -% Flag: 'c', 'b' oder 'cb' bzw. 'bc' +% Diese Funktion Zeichnet alle Vierecke mit ausgefuellten Flaechen. +% FLAG: +% c -> Koordinaten als rote Kreise darstellen +% b -> nur Kanten der Elemente einzeichnen +% n -> Normen auf ein Element einzeichnen (mit Laenge 1) +% a -> Normen auf ein Element einzeichnen (Laenge durch Flaecheninhalt) % % P.Schaefer @@ -25,7 +26,9 @@ elseif(optargin==1) if(ismember('c',varargin{1})) c = 1; end - if(ismember('n',varargin{1})) + if(ismember('a',varargin{1})) + n = 2; + elseif(ismember('n',varargin{1})) n = 1; end end @@ -55,7 +58,11 @@ else end if(n) - anorm = quad_norm(coordinates,elements); + if(n==2) + anorm = quad_norm(coordinates,elements,'w'); + else + anorm = quad_norm(coordinates,elements); + end for idx = 1:eles current = sum(coordinates(elements(idx,[2,3])',:),1)/2; current = [current ; current+anorm(idx,:)]; diff --git a/src/quad_norm.m b/src/quad_norm.m index 8fd0c78..6dbfe79 100644 --- a/src/quad_norm.m +++ b/src/quad_norm.m @@ -3,8 +3,9 @@ function n = quad_norm(coordinates, elements,varargin) % norm = tri_norm(coordinates, elements) % norm = tri_norm(coordinates, elements, 'FLAG') % -% This function builds the norm for every triangle -% +% Diese Funktion Berechnet die Orthogonalen mit Laenge 1 über alle Flächen +% FLAG: +% w -> Laenge entspricht Flaecheninhalt % P.Schaefer %% Parameterueberpruefung diff --git a/src/refineQuad.m b/src/refineQuad.m index 3c199fc..fffea50 100644 --- a/src/refineQuad.m +++ b/src/refineQuad.m @@ -11,24 +11,6 @@ for i = 1:c_loop if(type(i)==0) continue; elseif(type(i)==1) - newc = (coordinates(el(3),:)+coordinates(el(1),:))/2; - coordinates(c_coo+1,:) = newc; - elements(i,3) = c_coo+1; - - newc = newc+(coordinates(el(2),:)-coordinates(el(1),:)); - coordinates(c_coo+2,:) = newc; - - elements(c_ele+1,:) = [c_coo+1,c_coo+2,el(3)]; - elseif(type(i)==2) - newc = (coordinates(el(3),:)+coordinates(el(1),:))/2; - coordinates(c_coo+1,:) = newc; - elements(i,3) = c_coo+1; - - newc = newc+(coordinates(el(2),:)-coordinates(el(1),:)); - coordinates(c_coo+2,:) = newc; - - elements(c_ele+1,:) = [c_coo+1,c_coo+2,el(3)]; - elseif(type(i)==3) newc1 = (coordinates(el(2),:)+coordinates(el(1),:))/2; newc2 = (coordinates(el(3),:)+coordinates(el(1),:))/2; coordinates(c_coo+1,:) = newc1; @@ -45,6 +27,24 @@ for i = 1:c_loop coordinates(c_coo+4,:) = newc4; coordinates(c_coo+5,:) = newc5; elements(c_ele+3,:) = [c_coo+3,c_coo+4,c_coo+5]; + elseif(type(i)==2) + newc = (coordinates(el(3),:)+coordinates(el(1),:))/2; + coordinates(c_coo+1,:) = newc; + elements(i,3) = c_coo+1; + + newc = newc+(coordinates(el(2),:)-coordinates(el(1),:)); + coordinates(c_coo+2,:) = newc; + + elements(c_ele+1,:) = [c_coo+1,c_coo+2,el(3)]; + elseif(type(i)==3) + newc = (coordinates(el(3),:)+coordinates(el(1),:))/2; + coordinates(c_coo+1,:) = newc; + elements(i,3) = c_coo+1; + + newc = newc+(coordinates(el(2),:)-coordinates(el(1),:)); + coordinates(c_coo+2,:) = newc; + + elements(c_ele+1,:) = [c_coo+1,c_coo+2,el(3)]; end end end diff --git a/src/test_refine.m b/src/test_refine.m index 0446fa0..5fd7113 100644 --- a/src/test_refine.m +++ b/src/test_refine.m @@ -11,7 +11,7 @@ if(r_type > 3) marked = round(3*rand(1,eles)); - elseif(r_type < 0) + elseif(r_type <= 0) return; else marked = ones(1,eles)*r_type; -- 2.47.3