1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
|
! Copyright (C) 2006 Imperial College London and others.
!
! Please see the AUTHORS file in the main source directory for a full list
! of copyright holders.
!
! Prof. C Pain
! Applied Modelling and Computation Group
! Department of Earth Science and Engineering
! Imperial College London
!
! amcgsoftware@imperial.ac.uk
!
! This library is free software; you can redistribute it and/or
! modify it under the terms of the GNU Lesser General Public
! License as published by the Free Software Foundation; either
! version 2.1 of the License, or (at your option) any later version.
!
! This library is distributed in the hope that it will be useful,
! but WITHOUT ANY WARRANTY; without even the implied warranty of
! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
! Lesser General Public License for more details.
!
! You should have received a copy of the GNU Lesser General Public
! License along with this library; if not, write to the Free Software
! Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307
! USA
#include "fdebug.h"
module node_locking
use embed_python
use fields
use fldebug
use global_parameters, only : PYTHON_FUNC_LEN
use spud
implicit none
private
public :: get_locked_nodes
contains
subroutine get_locked_nodes(positions, locked_nodes, current_time)
!!< Return an array of nodes to be locked in mesh adaptivity. locked_nodes
!!< is allocated by this routine.
type(vector_field), intent(in) :: positions
integer, dimension(:), allocatable, intent(out) :: locked_nodes
real, optional, intent(in) :: current_time
character(len = *), parameter :: base_path = "/mesh_adaptivity/hr_adaptivity/node_locking"
character(len = PYTHON_FUNC_LEN) :: func
integer :: i, index, stat
integer, dimension(:), allocatable :: is_node_locked
real :: lcurrent_time
if(.not. have_option(base_path)) then
allocate(locked_nodes(0))
ewrite(2, *) "Number of locked nodes = 0"
return
end if
if(present(current_time)) then
lcurrent_time = current_time
else
call get_option("/timestepping/current_time", lcurrent_time, default = 0.0)
end if
call get_option(base_path // "/python", func)
allocate(is_node_locked(node_count(positions)))
call set_integer_array_from_python(func, len_trim(func), positions%dim, node_count(positions), &
& positions%val(1,:), positions%val(2,:), positions%val(3,:), lcurrent_time, &
& is_node_locked, stat)
if(stat /= 0) then
ewrite(-1, *) "Python error, Python string was:"
ewrite(-1, *) trim(func)
FLExit("Dying")
end if
allocate(locked_nodes(count(is_node_locked /= 0)))
ewrite(2, "(a,i0)") "Number of locked nodes = ", size(locked_nodes)
index = 0
do i = 1, size(is_node_locked)
if(is_node_locked(i) /= 0) then
index = index + 1
assert(index <= size(locked_nodes))
locked_nodes(index) = i
end if
end do
assert(index == size(locked_nodes))
deallocate(is_node_locked)
end subroutine get_locked_nodes
end module node_locking
|