From: wolf Date: Tue, 16 Feb 1999 17:15:02 +0000 (+0000) Subject: Add this script also. X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=3d659676bff7ea72d19eb78a3cdc9195dfac956b;p=dealii-svn.git Add this script also. git-svn-id: https://svn.dealii.org/trunk@818 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/source/fe/scripts/3d/lagrange-cubic b/deal.II/deal.II/source/fe/scripts/3d/lagrange-cubic new file mode 100644 index 0000000000..1716d22041 --- /dev/null +++ b/deal.II/deal.II/source/fe/scripts/3d/lagrange-cubic @@ -0,0 +1,161 @@ +# --------------------------------- For 3d --------------------------------- +# -- Use the following maple script to generate the basis functions, +# -- gradients and prolongation matrices as well as the mass matrix. +# -- Make sure that the files do not exists beforehand, since output +# -- is appended instead of overwriting previous contents. +# -- +# -- You should only have to change the very first lines for polynomials +# -- of higher order. +# +# -------------------------------------------------------------------------- +# +# $Id$ +# Author: Wolfgang Bangerth, 1999 + + read "lagrange-tools": + + n_functions := 64: + n_face_functions := 16: + n_constraints := 45: + + trial_function := ((a1 + a2*xi + a3*xi*xi + a4*xi*xi*xi) + + (b1 + b2*xi + b3*xi*xi + b4*xi*xi*xi)*eta + + (c1 + c2*xi + c3*xi*xi + c4*xi*xi*xi)*eta*eta + + (d1 + d2*xi + d3*xi*xi + d4*xi*xi*xi)*eta*eta*eta) + + ((e1 + e2*xi + e3*xi*xi + e4*xi*xi*xi) + + (f1 + f2*xi + f3*xi*xi + f4*xi*xi*xi)*eta + + (g1 + g2*xi + g3*xi*xi + g4*xi*xi*xi)*eta*eta + + (h1 + h2*xi + h3*xi*xi + h4*xi*xi*xi)*eta*eta*eta)*zeta + + ((i1 + i2*xi + i3*xi*xi + i4*xi*xi*xi) + + (j1 + j2*xi + j3*xi*xi + j4*xi*xi*xi)*eta + + (k1 + k2*xi + k3*xi*xi + k4*xi*xi*xi)*eta*eta + + (l1 + l2*xi + l3*xi*xi + l4*xi*xi*xi)*eta*eta*eta)*zeta*zeta + + ((m1 + m2*xi + m3*xi*xi + m4*xi*xi*xi) + + (n1 + n2*xi + n3*xi*xi + n4*xi*xi*xi)*eta + + (o1 + o2*xi + o3*xi*xi + o4*xi*xi*xi)*eta*eta + + (p1 + p2*xi + p3*xi*xi + p4*xi*xi*xi)*eta*eta*eta)*zeta*zeta*zeta: + face_trial_function := subs(zeta=0, trial_function): + # note: support_points[i] is a vector which is indexed from + # one and not from zero! + support_points := array(0..n_functions-1): + + support_points_fill_vertices (0, support_points): + support_points_fill_lines (8, 2, support_points): + support_points_fill_quads (32, 2, support_points): + support_points_fill_hex (56, 2, support_points): + + face_support_points := array(0..n_face_functions-1): + face_support_points[0] := [0,0]: + face_support_points[1] := [1,0]: + face_support_points[2] := [1,1]: + face_support_points[3] := [0,1]: + face_support_points[4] := [1/3,0]: + face_support_points[5] := [2/3,0]: + face_support_points[6] := [1,1/3]: + face_support_points[7] := [1,2/3]: + face_support_points[8] := [1/3,1]: + face_support_points[9] := [2/3,1]: + face_support_points[10] := [0,1/3]: + face_support_points[11] := [0,2/3]: + face_support_points[12] := [1/3,1/3]: + face_support_points[13] := [1/3,2/3]: + face_support_points[14] := [2/3,1/3]: + face_support_points[15] := [2/3,2/3]: + + # list of functions which are at face 0, used to compute + # the constraints on a face + constrained_face_function := array (0..n_face_functions-1): + # the list of points at which we want the functions at + # faces to be evaluated + constrained_face_support_points := array(0..n_constraints-1): + constrained_face_function[0] := 0: + constrained_face_function[1] := 1: + constrained_face_function[2] := 2: + constrained_face_function[3] := 3: + constrained_face_function[4] := 8: + constrained_face_function[5] := 9: + constrained_face_function[6] := 10: + constrained_face_function[7] := 11: + constrained_face_function[8] := 12: + constrained_face_function[9] := 13: + constrained_face_function[10] := 14: + constrained_face_function[11] := 15: + constrained_face_function[12] := 32: + constrained_face_function[13] := 33: + constrained_face_function[14] := 34: + constrained_face_function[15] := 35: + + constrained_face_support_points[0] := array(0..1, [1/2,1/2]): # center vertex + constrained_face_support_points[1] := array(0..1, [1/2,0]): # centers of large lines + constrained_face_support_points[2] := array(0..1, [1,1/2]): + constrained_face_support_points[3] := array(0..1, [1/2,1]): + constrained_face_support_points[4] := array(0..1, [0,1/2]): + constrained_face_support_points[5] := array(0..1, [1/2,1/6]): # lines from center to boundary + constrained_face_support_points[6] := array(0..1, [1/2,2/6]): + constrained_face_support_points[7] := array(0..1, [4/6,1/2]): + constrained_face_support_points[8] := array(0..1, [5/6,1/2]): + constrained_face_support_points[9] := array(0..1, [1/2,4/6]): + constrained_face_support_points[10] := array(0..1, [1/2,5/6]): + constrained_face_support_points[11] := array(0..1, [1/6,1/2]): + constrained_face_support_points[12] := array(0..1, [2/6,1/2]): + constrained_face_support_points[13] := array(0..1, [1/6,0]): # children of bounding lines + constrained_face_support_points[14] := array(0..1, [2/6,0]): + constrained_face_support_points[15] := array(0..1, [4/6,0]): + constrained_face_support_points[16] := array(0..1, [5/6,0]): + + constrained_face_support_points[17] := array(0..1, [1,1/6]): + constrained_face_support_points[18] := array(0..1, [1,2/6]): + constrained_face_support_points[19] := array(0..1, [1,4/6]): + constrained_face_support_points[20] := array(0..1, [1,5/6]): + + constrained_face_support_points[21] := array(0..1, [1/6,1]): + constrained_face_support_points[22] := array(0..1, [2/6,1]): + constrained_face_support_points[23] := array(0..1, [4/6,1]): + constrained_face_support_points[24] := array(0..1, [5/6,1]): + + constrained_face_support_points[25] := array(0..1, [0,1/6]): + constrained_face_support_points[26] := array(0..1, [0,2/6]): + constrained_face_support_points[27] := array(0..1, [0,4/6]): + constrained_face_support_points[28] := array(0..1, [0,5/6]): + + constrained_face_support_points[29] := array(0..1, [1/6,1/6]): # child quads + constrained_face_support_points[30] := array(0..1, [2/6,1/6]): + constrained_face_support_points[31] := array(0..1, [1/6,2/6]): + constrained_face_support_points[32] := array(0..1, [2/6,2/6]): + + constrained_face_support_points[33] := array(0..1, [4/6,1/6]): + constrained_face_support_points[34] := array(0..1, [5/6,1/6]): + constrained_face_support_points[35] := array(0..1, [4/6,2/6]): + constrained_face_support_points[36] := array(0..1, [5/6,2/6]): + + constrained_face_support_points[37] := array(0..1, [4/6,4/6]): + constrained_face_support_points[38] := array(0..1, [5/6,4/6]): + constrained_face_support_points[39] := array(0..1, [4/6,5/6]): + constrained_face_support_points[40] := array(0..1, [5/6,5/6]): + + constrained_face_support_points[41] := array(0..1, [1/6,4/6]): + constrained_face_support_points[42] := array(0..1, [2/6,4/6]): + constrained_face_support_points[43] := array(0..1, [1/6,5/6]): + constrained_face_support_points[44] := array(0..1, [2/6,5/6]): + + + # do the real work + read "lagrange": + + + + # write data to files + print ("writing data to files"): + readlib(C): + C(phi_polynom, filename=cubic3d_shape_value): + C(grad_phi_polynom, filename=cubic3d_shape_grad): + C(grad_grad_phi_polynom, filename=cubic3d_shape_grad_grad): + C(prolongation, filename=cubic3d_prolongation): + C(restriction, filename=cubic3d_restriction): + C(interface_constraints, filename=cubic3d_constraints): + C(real_space_points, optimized, filename=cubic3d_real_points): + + writeto (cubic3d_unit_support_points): + print (support_points): + + \ No newline at end of file