Mercurial > repos > public > sbplib
changeset 588:42124009f940 feature/better_multiblock_defs
Add domain definition for circle.
author | Jonatan Werpers <jonatan@werpers.com> |
---|---|
date | Mon, 11 Sep 2017 13:50:29 +0200 |
parents | 25fdc7a625b6 |
children | 53d422efaabd |
files | +multiblock/+domain/Circle.m |
diffstat | 1 files changed, 98 insertions(+), 0 deletions(-) [+] |
line wrap: on
line diff
diff -r 25fdc7a625b6 -r 42124009f940 +multiblock/+domain/Circle.m --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/+multiblock/+domain/Circle.m Mon Sep 11 13:50:29 2017 +0200 @@ -0,0 +1,98 @@ +classdef Circle < multiblock.DefCurvilinear + properties + r, c + + hs + r_arc + omega + end + + methods + function obj = Circle(r, c, hs) + default_arg('r', 1); + default_arg('c', [0; 0]); + default_arg('hs', 0.435); + + + % alpha = 0.75; + % hs = alpha*r/sqrt(2); + + % Square should not be a square, it should be an arc. The arc radius + % is chosen so that the three angles of the meshes are all equal. + % This gives that the (half)arc opening angle of should be omega = pi/12 + omega = pi/12; + r_arc = hs*(2*sqrt(2))/(sqrt(3)-1); % = hs* 1/sin(omega) + c_arc = c - [(1/(2-sqrt(3))-1)*hs; 0]; + + cir = parametrization.Curve.circle(c,r,[-pi/4 pi/4]); + + c2 = cir(0); + c3 = cir(1); + + s1 = [-hs; -hs]; + s2 = [ hs; -hs]; + s3 = [ hs; hs]; + s4 = [-hs; hs]; + + sp2 = parametrization.Curve.line(s2,c2); + sp3 = parametrization.Curve.line(s3,c3); + + Se1 = parametrization.Curve.circle(c_arc,r_arc,[-omega, omega]); + Se2 = Se1.rotate(c,pi/2); + Se3 = Se2.rotate(c,pi/2); + Se4 = Se3.rotate(c,pi/2); + + + S = parametrization.Ti(Se1,Se2,Se3,Se4).rotate_edges(-1); + + A = parametrization.Ti(sp2, cir, sp3.reverse, Se1.reverse); + B = A.rotate(c,1*pi/2).rotate_edges(-1); + C = A.rotate(c,2*pi/2).rotate_edges(-1); + D = A.rotate(c,3*pi/2).rotate_edges(0); + + blocks = {S,A,B,C,D}; + blocksNames = {'S','A','B','C','D'}; + + conn = cell(5,5); + conn{1,2} = {'e','w'}; + conn{1,3} = {'n','s'}; + conn{1,4} = {'w','s'}; + conn{1,5} = {'s','w'}; + + conn{2,3} = {'n','e'}; + conn{3,4} = {'w','e'}; + conn{4,5} = {'w','s'}; + conn{5,2} = {'n','s'}; + + boundaryGroups = struct(); + boundaryGroups.E = multiblock.BoundaryGroup({2,'e'}); + boundaryGroups.N = multiblock.BoundaryGroup({3,'n'}); + boundaryGroups.W = multiblock.BoundaryGroup({4,'n'}); + boundaryGroups.S = multiblock.BoundaryGroup({5,'e'}); + boundaryGroups.all = multiblock.BoundaryGroup({{2,'e'},{3,'n'},{4,'n'},{5,'e'}}); + + obj = obj@multiblock.DefCurvilinear(blocks, conn, boundaryGroups, blocksNames); + + obj.r = r; + obj.c = c; + obj.hs = hs; + obj.r_arc = r_arc; + obj.omega = omega; + end + + function ms = getGridSizes(obj, m) + m_S = m; + + % m_Radial + s = 2*obj.hs; + innerArc = obj.r_arc*obj.omega; + outerArc = obj.r*pi/2; + shortSpoke = obj.r-s/sqrt(2); + x = (1/(2-sqrt(3))-1)*obj.hs; + longSpoke = (obj.r+x)-obj.r_arc; + m_R = parametrization.equal_step_size((innerArc+outerArc)/2, m_S, (shortSpoke+longSpoke)/2); + + ms = {[m_S m_S], [m_R m_S], [m_S m_R], [m_S m_R], [m_R m_S]}; + end + end +end