• R/O
  • SSH

treelm: 提交

Repository of the treelm library. Now found at https://github.com/apes-suite/treelm


Commit MetaInfo

修訂c2eaab17a4e20504605d5af9b595010ba24918f3 (tree)
時間2023-05-09 18:20:32
作者Kannan Masilamani <kannan.masilamani@dlr....>
CommiterKannan Masilamani

Log Message

Moved search for neighbor when a point is outside domain from
tem_cano_storePntsInSubTree to tem_cano_initSubTree.
In tem_cano_storePntsInSubTree, a treeID of a point is searched in subTreeIds
so a point outside domain should be found since its neighbor is added to subTree
map2global.

Change Summary

差異

diff -r be8e137a8a01 -r c2eaab17a4e2 source/tem_canonical_module.f90
--- a/source/tem_canonical_module.f90 Mon May 08 17:11:12 2023 +0200
+++ b/source/tem_canonical_module.f90 Tue May 09 11:20:32 2023 +0200
@@ -722,8 +722,9 @@
722722 !> If true then elements not intersected are added to subTree
723723 logical, intent(in) :: shapeInverted
724724 ! --------------------------------------------------------------------------
725- integer :: iCano, tLevel, elemPos, iElem, dPos, maxLevel
726- integer(kind=long_k) :: treeID
725+ integer :: iCano, tLevel, elemPos, iElem, dPos, maxLevel, iQQN
726+ integer :: neighPos, coordOfId(4)
727+ integer(kind=long_k) :: treeID, tOffset, neighID
727728 logical :: wasAdded, intersects, addToSubTree
728729 type(tem_cube_type) :: cube
729730 ! --------------------------------------------------------------------------
@@ -756,6 +757,31 @@
756757 tLevel = tem_levelOf( treeID )
757758 countElems( tLevel ) = countElems( tLevel ) + 1
758759 end if ! wasAdded
760+ else
761+ ! add a neighbor element of this point if any of the neighbor exist to
762+ ! interpolate to a point
763+ coordOfId = tem_CoordOfId( treeID )
764+ tOffset = tem_FirstIdAtLevel( coordOfId(4) )
765+ directionLoop: do iQQN = 1, qQQQ
766+ neighID = tem_IdOfCoord( &
767+ & [ coordOfId(1) + qOffset( iQQN, 1 ), &
768+ & coordOfId(2) + qOffset( iQQN, 2 ), &
769+ & coordOfId(3) + qOffset( iQQN, 3 ), &
770+ & coordOfId(4) ], tOffset)
771+ neighPos = tem_PosOfId( neighID, inTree%treeID)
772+ if (neighPos > 0) then
773+ call append( me = map2global, &
774+ & pos = dpos, &
775+ & val = neighPos, &
776+ & wasAdded = wasAdded )
777+
778+ ! Count up if it was added
779+ if (wasAdded) then
780+ tLevel = tem_levelOf( neighID )
781+ countElems( tLevel ) = countElems( tLevel ) + 1
782+ end if ! wasAdded
783+ end if ! neighPos > 0
784+ end do directionLoop
759785 end if ! elemPos > 0
760786 else
761787 do iElem = 1, inTree%nElems
@@ -818,11 +844,11 @@
818844 !> growing array to store tracking points
819845 type(tem_grwPoints_type), intent(inout) :: grwPnts
820846 ! --------------------------------------------------------------------------
821- integer :: nElems, nPoints, maxLevel, elemPos, neighPos, coordOfId(4)
822- integer :: iCano, iPnt, iQQN
847+ integer :: nElems, nPoints, maxLevel, elemPos
848+ integer :: iCano, iPnt
823849 real(kind=rk) :: coord(3), offset_a, offset_b, offset_c
824850 real(kind=rk) :: unit_vec_a(3),unit_vec_b(3), unit_vec_c(3)
825- integer(kind=long_k) :: treeID, tOffset, neighID
851+ integer(kind=long_k) :: treeID
826852 integer(kind=long_k), allocatable :: subTreeID(:)
827853 ! --------------------------------------------------------------------------
828854 maxLevel = inTree%global%maxLevel
@@ -866,28 +892,6 @@
866892 & val = coord )
867893
868894 countPoints = countPoints + 1
869- else
870- ! add this if any of the neighbor exist to interpolate/extrapolate
871- coordOfId = tem_CoordOfId( treeID )
872- tOffset = tem_FirstIdAtLevel( coordOfId(4) )
873- directionLoop: do iQQN = 1, qQQQ
874- neighID = tem_IdOfCoord( &
875- & [ coordOfId(1) + qOffset( iQQN, 1 ), &
876- & coordOfId(2) + qOffset( iQQN, 2 ), &
877- & coordOfId(3) + qOffset( iQQN, 3 ), &
878- & coordOfId(4) ], tOffset)
879- neighPos = tem_PosOfId( neighID, inTree%treeID)
880- if (neighPos > 0) then
881- ! append the physical points to the growing array of points
882- call append( me = grwpnts, &
883- & val = coord )
884-
885- countPoints = countPoints + 1
886-
887- exit directionLoop
888- end if
889- end do directionLoop
890-
891895 end if
892896 end do !iPoint
893897 end do !iCano
diff -r be8e137a8a01 -r c2eaab17a4e2 source/tem_shape_module.f90
--- a/source/tem_shape_module.f90 Mon May 08 17:11:12 2023 +0200
+++ b/source/tem_shape_module.f90 Tue May 09 11:20:32 2023 +0200
@@ -548,7 +548,7 @@
548548
549549 countPoints = 0
550550 call init(me=grwPnts)
551- if ( storePnts ) then
551+ if ( storePnts .and. (map2global%nVals > 0) ) then
552552 select case (trim(me%kind))
553553 case ('canoND')
554554 call tem_cano_storePntsInSubTree( me = me%canoND(:), &
Show on old repository browser