This subroutine calculates the surface area attached to each point
Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
type(tem_surfData_type), | intent(inout) | :: | me |
datatype to store the surface information |
subroutine tem_calcTriaAreas( me ) ! --------------------------------------------------------------------------- !> datatype to store the surface information type( tem_surfData_type ), intent(inout) :: me ! --------------------------------------------------------------------------- ! counters integer :: iTria ! temporary variable to store the triangle area real(kind=rk) :: area ! temporary vectors of two sides of the triangle real(kind=rk) :: vec1(3), vec2(3) ! temporary variable storing the start and end positions in the pointCoords ! array integer :: minPos1, maxPos1, minPos2, maxPos2, minPos3, maxPos3 ! --------------------------------------------------------------------------- ! reset the array of surface areas me%surfArea = 0._rk ! loop over the triangles and ... do iTria = 1, me%nTrias ! ... set the minimum and maximum positions in the pointCoords array minPos1 = (me%trias(1, iTria) - 1 )*3 + 1 maxPos1 = (me%trias(1, iTria) - 1 )*3 + 3 minPos2 = (me%trias(2, iTria) - 1 )*3 + 1 maxPos2 = (me%trias(2, iTria) - 1 )*3 + 3 minPos3 = (me%trias(3, iTria) - 1 )*3 + 1 maxPos3 = (me%trias(3, iTria) - 1 )*3 + 3 ! if( me%parentIDs(iLevel)%ptrs(me%trias(1, iTria)) > 0 .and. & ! & me%parentIDs(iLevel)%ptrs(me%trias(2, iTria)) > 0 .and. & ! & me%parentIDs(iLevel)%ptrs(me%trias(3, iTria)) > 0 )then ! ... calculate two vectors of the triangle vec1 = me%pointCoords(minPos1:maxPos1) & & - me%pointCoords(minPos2:maxPos2) vec2 = me%pointCoords(minPos1:maxPos1) & & - me%pointCoords(minPos3:maxPos3) ! ... calculate the surface area of the triangle area = sqrt( sum( cross_product3D( vec1, vec2 )**2)) ! ... distribute it among the attached surface points me%surfArea( me%trias(1, iTria)) = me%surfArea( me%trias(1, iTria)) & & + area/3._rk me%surfArea( me%trias(2, iTria)) = me%surfArea( me%trias(2, iTria)) & & + area/3._rk me%surfArea( me%trias(3, iTria)) = me%surfArea( me%trias(3, iTria)) & & + area/3._rk ! end if end do end subroutine tem_calcTriaAreas