restart: with(plots): with(plottools): # A global example: # # A fluid flows within a cylindrical duct, passes through a diaphragm drilled by several hole, # enters a diffuser of complex shape wich raccords to a rectangular bow wich contains # heat exchanger made of 4 bundles of 10 plates. # # This is part of a real problem whose goal is to find a topology of the diaphragm (positions # and diameters of the whole) which minimizes the differences in heat flow within each # pairs of plates). # This is achieved (still under development) with the GlobalOptimization Toolbox and (for the moment) # "zero order" hydraulic models. The real problem is based on CFD simulations made by # a commercial code. # # The purpose is to demonstrate here the many capabilities of Maple in drawing realistic complex # systems. # # Some variable names being french written here are a few translations : # conduite : duct # diffuseur : diffuser # diaphragme : diaphragm # faisceau : bundle # tourner : (to) rotate # complement : complement # dilate : (to) dilate # echec : fiailure (fail) # # Duct P := 8: Q := P+1; R := P+1+1; conduite := cylinder([0, 0, -20], Q, 20, color=gray, capped=false, strips=200, style=surface, transparency=0.4): a := display(annulus(Q..Q+0.01), scaling=constrained): conduite := extrude(a, -20..0, color=gray, transparency=0.4): # The "morphing" function f is aimed to continuously transform a circle into a square. # It is base on the fact that the graph of the relation x^n+y^n-1 =0 is a circle # for n=2 and tends to the unit square as n tends to infinity f := u -> 2+9*(u/10)^3: diffuseur := implicitplot3d((abs(x))^f(z)+(abs(y))^f(z)-Q^f(z), x=-Q..Q, y=-Q..Q, z=0..10, style=surface, color=red, grid=[10, 100, 100], scaling=constrained, transparency=0.4): # Heat exchanger bloc2D := implicitplot((abs(x))^f(10)+(abs(y))^f(10)-Q^f(10), x=-Q..Q, y=-Q..Q, grid=[100, 100], color=gold): bloc3D := display(extrude(bloc2D, 10..30), transparency=0.4): # Bundles K := 10: faisceau_1 := display( seq( polygon( [ [4/P+(k-1)*P/K, 1/2, 12], [4/P+(k-1)*P/K, P , 12], [4/P+(k-1)*P/K, P , 30], [4/P+(k-1)*P/K, 1/2, 30] ], color=gray ), k=1..K+1 ), scaling=constrained, axes=none ): faisceau := display( seq( polygon( [ [4/P+(k-1)*P/K, 1/2, 12], [4/P+(k-1)*P/K, P , 12], [4/P+(k-1)*P/K, P , 30], [4/P+(k-1)*P/K, 1/2, 30] ], color=gray ), k=1..K+1 ), scaling=constrained, axes=none ): tourner := a -> transform((x, y, z) -> [x*cos(a)+y*sin(a), -x*sin(a)+y*cos(a), z]): IiIq IiM1 # Diaphragm # # A two-step construction of the diaphragm: # 1/ The pattern (a square drilled with a circular hole) is first constructed (name "a"). # 2/ This pattern is periodically replicated to fill (as much as possible) the unit disc (name "diaphragme"). m := 0.2: eps := 0.001: ij := [ seq( seq( `if`( (i+m/2)^2+(j+m/2)^2 < 1 and (i+m/2)^2+(j-m/2)^2 < 1 and (i-m/2)^2+(j-m/2)^2 < 1 and (i-m/2)^2+(j+m/2)^2 < 1, [i,j], NULL), i in [seq(-1..1, m)] ), j in [seq(-1..1, m)] ) ] : a := polygon((m/2) *~ [ [-1, -1], [1, -1], [1, -eps], [1/2, -eps], seq([1/2*cos(-t), 1/2*sin(-t)], t in [seq(eps..2*Pi-eps, (2*Pi-2*eps)/100)]), [1/2, +eps], [1,+eps], [1, 1], [-1, 1], [-1, -1] ], style='polygon' ): T := (X,Y) -> transform((x, y) -> [x+X, y+Y]): b := display(seq(T(op(u))(a), u in ij)): diaphragme := extrude(b, 0..0.2, color=blue, transparency=0): display( diaphragme, scaling=constrained ): # Extension of the diaphragm to the unit circle # # The complementary polygon is termed "complement" IJ := copy(ij): Go_North := proc(IJ, P) local S, lop, pol: lop := map(u -> if u=P and u > P then u end if, IJ); if numelems(lop) >=1 then pol := sort(lop , key=(x->x)); else pol := NULL: end if: end proc: Go_East := proc(IJ, P) local S, lop, pol: lop := map(u -> if u=P and u > P then u end if, IJ); if numelems(lop) >=1 then pol := sort(lop , key=(x->x)); else pol := NULL: end if: end proc: Go_South := proc(IJ, P) local S, lop, pol: lop := map(u -> if u=P and u < P then u end if, IJ); if numelems(lop) >=1 then pol := sort(lop , key=(x->x))[-1]; else pol := NULL: end if: end proc: Go_West := proc(IJ, P) local S, lop, pol: lop := map(u -> if u=P and u < P then u end if, IJ); if numelems(lop) >=1 then pol := sort(lop , key=(x->x))[-1]; else pol := NULL: end if: end proc: __GO := proc(__to, IJ, P) #DEBUG(): if __to = "north" then Go_North(IJ, P) elif __to = "east" then Go_East (IJ, P) elif __to = "south" then Go_South(IJ, P) elif __to = "west" then Go_West (IJ, P) end if; end proc: min1 := min(op~(IJ)): POL := [ sort( map(u -> if u = min1 then u end if, IJ) ) ]: dir := ["north", "east", "south", "west", "north"]: first := 1: new := NULL: while new <> POL do echec := false: new := __GO(dir[first], IJ, POL[-1]); if new <> NULL then POL := [POL[], new]; else echec := true: new := __GO(dir[first+1], IJ, POL[-1]); if new <> NULL then POL := [POL[], new]; else if echec then first := first+1: else echec := false: end if: end if: end if: end do: POL: pol := NULL: for p in POL do pol := pol, p +~ (m/2) * signum~(p) end do: pol := [pol]: max1 := max(op~(1, pol)): first := sort( map(u -> if u = max1 then u end if, pol) )[-1]: pos := ListTools:-Search(first, pol): rpol := ListTools:-Rotate(pol, pos-1): theta := arctan(first/first): eps := 0.001: complement := [ seq([cos(t+theta), sin(t+theta)], t in [seq(eps..2*Pi-eps, (2*Pi-2*eps)/100)]), rpol[], rpol, [cos(eps+theta), sin(eps+theta)] ]: complement := display(polygon(complement, style='polygon')): complement := extrude(complement, 0..0.2, color=blue, transparency=0): # Assembly of the different elements. # # The image is not display for the file would have been too big (28 Mo) dilate := transform((x, y, z) -> [Q*x, Q*y, 10*z-4]): display( conduite, dilate(complement), dilate(diaphragme), diffuseur, bloc3D, faisceau, tourner(Pi/2)(faisceau), tourner(Pi)(faisceau), tourner(3*Pi/2)(faisceau) ): Here are a few extra constructions that could be useful to anyone interested # "Conical" extrusion of an annulus a := display(annulus(8..8.5)): b := extrude(a, 0..40, color=gray, transparency=0.4): T := plottools:-transform((x, y, z) -> [x*(1+z/40), y*(1+z/40), z]): display(T(b), scaling=constrained): # Different extrusions of an annulus morphed into another shape a := display(annulus(8..8.5)): b := extrude(a, 0..40, color=gray, transparency=0.): Tf := u -> 1+9*(u/40)^2: T := plottools:-transform((x, y, z) -> [signum(x)*(abs(x)/8)^Tf(z), signum(y)*(abs(y)/8)^Tf(z), z]): display(T(b), scaling=unconstrained): a := display(annulus(8..8.5)): b := extrude(a, 0..40, color=gray, transparency=0.): Tf := u -> 1+2*(u/40): T := plottools:-transform((x, y, z) -> [signum(y)*(abs(x)/8.5)^Tf(z), signum(x)*(abs(y)/8.5)^Tf(z), z]): display(T(b), scaling=unconstrained): # From circle to square (from S1 in L2 norm to S1 in L^16 norm) a := display(annulus(8..8.5)): b := extrude(a, 0..40, color=gray, transparency=0.): r0 := 8.25: Tf := u -> 1/(1+(u/20)^4): T := plottools:-transform((x, y, z) -> [signum(x)*(abs(x)/r0)^Tf(z), signum(y)*(abs(y)/r0)^Tf(z), z]): display(T(b), grid=[100\$3], scaling=unconstrained): # Extrusion of a "squared annulus" eps := 0.001: a := display( polygon( [ [-1, -1], [1, -1], [1, 1], [-1, 1], [-1, eps], [-1/2, eps], [-1/2, 1/2], [1/2, 1/2], [1/2, -1/2], [-1/2, -1/2], [-1/2, -eps], [-1, -eps], [-1, -1] ], style='polygon' ) ): b := extrude(a, 0..2, color=gray, transparency=0): display(b, scaling=constrained): # Box with a cylindrical hole eps := 0.001: a := (x,y) -> polygon( [ [-1+x, -1+y], [1+x, -1+y], [1+x, -eps+y], [1/2+x, -eps+y], seq([1/2*cos(-t)+x, 1/2*sin(-t)+y], t in [seq(eps..2*Pi-eps, (2*Pi-2*eps)/100)]), [1/2+x, +eps+y], [1+x,+eps+y], [1+x, 1+y], [-1+x, 1+y], [-1+x, -1+y] ], style='polygon' ): b := extrude(display(a(0, 0)), 0..0.5, color=gray, transparency=0): display(b, scaling=constrained): # Extrusion of a grid of "circularly holed squared annuli" b := display(seq(seq(a(i,j), i in [seq(-4..4, 2)]), j in [seq(-4..4, 2)])): diaphragme := extrude(b, -10..-8, color=blue, transparency=0): display( diaphragme, scaling=constrained ):