1 | include(dom.inc)
2 | SUBROUTINE calculatearea(nb_nodes, nodelist, area)
3 |
4 | ! ================================================================!
5 | ! !
6 | ! calculatearea.F : Calculates the surface of any kind of face !
7 | ! !
8 | ! in : The number of nodes of the face 'nb_nodes' !
9 | ! The id list of the nodes !
10 | ! out : The face's surface value 'area' !
11 | ! !
12 | ! author : J. AMAYA (avril 2007) !
13 | ! !
14 | ! ================================================================!
15 |
16 | USE datas
17 | IMPLICIT NONE
18 |
19 | ! IN
20 | DOM_INT :: nb_nodes
21 | DOM_INT, DIMENSION(nb_nodes) :: nodelist
22 |
23 | ! OUT
24 | DOM_REAL :: area
25 |
26 | ! LOCAL
27 | DOM_REAL :: s, a, b, c, a1, a2
28 | DOM_REAL :: x1, y1, z1, x2, y2, z2, x3, y3, z3
29 | DOM_INT :: ierr
30 |
31 | ! -------------------!
32 | ! For the edge in 2D !
33 | ! -------------------!
34 | IF (nb_nodes.eq.2) THEN
35 | x1 = node_list(1, nodelist(1))
36 | y1 = node_list(2, nodelist(1))
37 |
38 | x2 = node_list(1, nodelist(2))
39 | y2 = node_list(2, nodelist(2))
40 |
41 | a = x2 - x1
42 | b = y2 - y1
43 |
44 | area = SQRT( a*a + b*b )
45 | ! -------------------!
46 | ! For the triangle : !
47 | ! -------------------!
48 | ELSE IF (nb_nodes.eq.3) THEN
49 |
50 | x1 = node_list(1, nodelist(1))
51 | y1 = node_list(2, nodelist(1))
52 | z1 = node_list(3, nodelist(1))
53 |
54 | x2 = node_list(1, nodelist(2))
55 | y2 = node_list(2, nodelist(2))
56 | z2 = node_list(3, nodelist(2))
57 |
58 | x3 = node_list(1, nodelist(3))
59 | y3 = node_list(2, nodelist(3))
60 | z3 = node_list(3, nodelist(3))
61 |
62 | a=SQRT((x2-x1)*(x2-x1)+(y2-y1)*(y2-y1)+(z2-z1)*(z2-z1))
63 | b=SQRT((x3-x2)*(x3-x2)+(y3-y2)*(y3-y2)+(z3-z2)*(z3-z2))
64 | c=SQRT((x1-x3)*(x1-x3)+(y1-y3)*(y1-y3)+(z1-z3)*(z1-z3))
65 |
66 | s = (a + b + c) / 2
67 | area = SQRT(s*(s-a)*(s-b)*(s-c))
68 |
69 | ! ----------------!
70 | ! For the square: !
71 | ! ----------------!
72 | ELSE
73 |
74 | x1 = node_list(1, nodelist(1))
75 | y1 = node_list(2, nodelist(1))
76 | z1 = node_list(3, nodelist(1))
77 |
78 | x2 = node_list(1, nodelist(2))
79 | y2 = node_list(2, nodelist(2))
80 | z2 = node_list(3, nodelist(2))
81 |
82 | x3 = node_list(1, nodelist(3))
83 | y3 = node_list(2, nodelist(3))
84 | z3 = node_list(3, nodelist(3))
85 |
86 | a=SQRT((x2-x1)*(x2-x1)+(y2-y1)*(y2-y1)+(z2-z1)*(z2-z1))
87 | b=SQRT((x3-x2)*(x3-x2)+(y3-y2)*(y3-y2)+(z3-z2)*(z3-z2))
88 | c=SQRT((x1-x3)*(x1-x3)+(y1-y3)*(y1-y3)+(z1-z3)*(z1-z3))
89 |
90 | s = (a + b + c) / 2
91 | a1 = SQRT(s*(s-a)*(s-b)*(s-c))
92 |
93 | x2 = node_list(1, nodelist(4))
94 | y2 = node_list(2, nodelist(4))
95 | z2 = node_list(3, nodelist(4))
96 |
97 | a=SQRT((x2-x1)*(x2-x1)+(y2-y1)*(y2-y1)+(z2-z1)*(z2-z1))
98 | b=SQRT((x3-x2)*(x3-x2)+(y3-y2)*(y3-y2)+(z3-z2)*(z3-z2))
99 |
100 | s = (a + b + c) / 2
101 | a2 = SQRT(s*(s-a)*(s-b)*(s-c))
102 |
103 | area = a1 + a2
104 |
105 | ENDIF
106 |
107 | END SUBROUTINE calculatearea
calculatearea.F could be called by: