FrontISTR 5.2.0
Large-scale structural analysis program with finit element method
Loading...
Searching...
No Matches
ttable.f90
Go to the documentation of this file.
1!-------------------------------------------------------------------------------
2! Copyright (c) 2019 FrontISTR Commons
3! This software is released under the MIT License, see LICENSE.txt
4!-------------------------------------------------------------------------------
7module m_table
8 implicit none
9 integer, parameter, private :: kreal = kind(0.0d0)
10
11 !
12 ! The length of the keys
13 !
14 integer, parameter :: dict_key_length = 128
15 integer, parameter :: maxindex = 20
16 integer, parameter :: maxnval = 10000
17
18 !
19 ! The data that will be stored with each key
20 type ttable
21 integer :: ndepends
22 integer :: tbcol, tbrow
23 integer :: tbindex(maxindex)
24 real(kind=kreal), pointer :: tbval(:,:)=>null()
25 end type ttable
26
27 !
28 ! The "null" value for these data
29 !
30 type(ttable), parameter :: dict_null = ttable( 0,0,0, 0,null() )
31
32 interface assignment(=)
33 module procedure tablecopy
34 end interface
35
36 interface operator(==)
37 module procedure tablecompare
38 end interface
39
40contains
41
42 subroutine init_table( table, ndp, col, row, tbval )
43 type( ttable ), intent(inout) :: table
44 integer, intent(in) :: ndp, col, row
45 real(kind=kreal), intent(in) :: tbval(col,row)
46 integer :: i,j
47 real(kind=kreal) :: tbindexval(maxnval)
48 if( associated( table%tbval ) ) deallocate( table%tbval )
49 table%tbcol = col
50 table%tbrow = row
51 table%ndepends = ndp
52 allocate( table%tbval( col, row ) )
53 do i=1, col
54 do j=1, row
55 table%tbval(i,j) = tbval(i,j)
56 enddo
57 enddo
58
59 do i=1,ndp
60 table%tbindex(i) =1
61 tbindexval(1) = table%tbval(table%tbcol-i+1, 1)
62 ! print *, 0, i, table%tbindex(i), tbindexval(1)
63
64 do j=2, table%tbrow
65 if( any( tbindexval(1:table%tbindex(i))==table%tbval(table%tbcol-i+1, j) ) ) cycle
66 table%tbindex(i) = table%tbindex(i)+1
67 tbindexval(table%tbindex(i)) = table%tbval(table%tbcol-i+1, j)
68 ! print *, 1, i, table%tbindex(i), tbindexval(table%tbindex(i))
69 enddo
70 enddo
71
72 j=1
73 do i=1,ndp
74 j=j*table%tbindex(i)
75 enddo
76 if( j/= row) stop "Error in table defnition!"
77 ! print *, j,row, table%tbindex(1:ndp); pause
78 end subroutine
79
80 subroutine finalize_table( table )
81 type( ttable ), intent(inout) :: table
82 if( associated( table%tbval ) ) deallocate( table%tbval )
83 end subroutine
84
85 subroutine print_table( table, fname )
86 type( ttable ), intent(in) :: table
87 integer, intent(in) :: fname
88 integer :: i,j
89 write(fname,*) table%ndepends, table%tbcol, table%tbrow
90 do i=1,table%tbrow
91 write(fname,*) i,(table%tbval(j,i),j=1,table%tbcol)
92 enddo
93 end subroutine
94
95 logical function tablecompare( lhs, rhs )
96 type(ttable), intent(in) :: lhs
97 type(ttable), intent(in) :: rhs
98 integer :: i,j
99 tablecompare = .false.
100 if( lhs%ndepends /= rhs%ndepends ) return
101 if( lhs%tbcol /= rhs%tbcol ) return
102 if( lhs%tbrow /= rhs%tbrow ) return
103 do i=1, rhs%ndepends
104 if( lhs%tbindex(i) /= rhs%tbindex(i) ) return
105 enddo
106 do i=1, rhs%tbcol
107 do j=1, rhs%tbrow
108 if( lhs%tbval(i,j) /= rhs%tbval(i,j) ) return
109 enddo
110 enddo
111 tablecompare = .true.
112 end function
113
114 subroutine tablecopy( lhs, rhs )
115 type(ttable), intent(out) :: lhs
116 type(ttable), intent(in) :: rhs
117
118 lhs%ndepends = rhs%ndepends
119 lhs%tbcol = rhs%tbcol
120 lhs%tbrow = rhs%tbrow
121 lhs%tbindex(:) = rhs%tbindex(:)
122 if( associated( lhs%tbval ) ) deallocate( lhs%tbval )
123 if( rhs%tbcol<=0 .or. rhs%tbrow<=0 ) return
124 allocate( lhs%tbval( lhs%tbcol, lhs%tbrow ) )
125 lhs%tbval(:,:) = rhs%tbval(:,:)
126 end subroutine
127
128end module m_table
129
130!======================================================================!
131!
137!
138!======================================================================!
139
141 use m_table, dict_data => ttable
142 implicit none
143 integer, parameter, private :: kreal = kind(0.0d0)
144
145 private :: gettablegrad, gettabledata
146
147 include "dictionary.f90"
148
149
152 subroutine fetch_table( key, dict, dicval, ierr )
153 character(len=*), intent(in) :: key
154 type(dict_struct), pointer :: dict
155 logical, intent(out) :: ierr
156 type(dict_data), pointer :: dicval
157
158 dicval => dict_get_key( dict, key )
159 ierr = .false.
160 if( .not. associated(dicval) ) then
161 ierr=.true.; return
162 endif
163
164 end subroutine
165
166
168 integer function fetch_tablerow( key, dict )
169 character(len=*), intent(in) :: key
170 type(dict_struct), pointer :: dict
171 type(dict_data), pointer :: dicval
172
173 dicval => dict_get_key( dict, key )
174 fetch_tablerow = -1
175 if( .not. associated(dicval) ) return
176 fetch_tablerow = dicval%tbrow
177 ! call finalize_table( dicval )
178 end function
179
182 subroutine fetch_tablegrad( key, a, dict, outa, ierr )
183 character(len=*), intent(in) :: key
184 real(kind=kreal), intent(in) :: a(:)
185 type(dict_struct), pointer :: dict
186 real(kind=kreal), intent(out) :: outa
187 logical, intent(out) :: ierr
188
189 type(dict_data), pointer :: dicval
190 integer :: na, dd, crow, cindex
191 dicval => dict_get_key( dict, key )
192 ierr = .false.
193 if( .not. associated(dicval) ) then
194 ierr=.true.; return
195 endif
196
197 cindex = 1; crow=1
198 dd = dicval%tbrow
199 na = 1
200 if( size(a) > dicval%ndepends ) then ! in case no dependent definition
201 na = size(a) - dicval%ndepends+1
202 endif
203
204 call gettablegrad( a(na:), cindex, dicval, dd, crow, outa )
205 ! call finalize_table( dicval )
206
207 end subroutine
208
210 recursive subroutine gettablegrad( a, cindex, table, dd, crow, outa )
211 real(kind=kreal), intent(in) :: a(:)
212 integer, intent(inout) :: cindex
213 type(dict_data) :: table
214 integer, intent(inout) :: dd, crow
215 real(kind=kreal), intent(out) :: outa
216
217 integer :: i, ccol, ddd
218 real(kind=kreal) :: val1, val2, lambda
219
220 ddd = dd / table%tbindex(cindex)
221 ccol = table%tbcol-cindex+1
222
223 if( ddd==1 ) then
224 if( a(cindex)<table%tbval(2, crow) ) then
225 outa = 0.d0
226 elseif( a(cindex)>=table%tbval(2, crow+dd-1) ) then
227 outa = 0.d0
228 else
229 do i=crow, crow+dd-2
230 if( a(cindex)>=table%tbval(2, i) .and. a(cindex)<table%tbval(2, i+1) ) then
231 outa = (table%tbval(1, i+1)-table%tbval(1, i))/(table%tbval(2, i+1)-table%tbval(2, i))
232 exit
233 endif
234 enddo
235 endif
236 else
237 if( a(cindex)<=table%tbval(ccol, crow) ) then
238 cindex = cindex+1
239 dd = ddd
240 call gettablegrad( a, cindex, table, dd, crow, outa )
241 elseif( a(cindex)>=table%tbval(ccol, crow+dd-1) ) then
242 crow = crow+dd-ddd
243 cindex = cindex+1
244 dd = ddd
245 call gettablegrad( a, cindex, table, dd, crow, outa )
246 else
247 do i=crow, crow+dd-2, ddd
248 if( a(cindex)==table%tbval(ccol, i) ) then
249 crow = i
250 cindex = cindex+1
251 dd = ddd
252 call gettablegrad( a, cindex, table, dd, crow, outa )
253 exit
254 elseif( a(cindex)==table%tbval(ccol, i+ddd) ) then
255 crow = i+ddd
256 cindex = cindex+1
257 dd = ddd
258 call gettablegrad( a, cindex, table, dd, crow, outa )
259 exit
260 elseif( a(cindex)>table%tbval(ccol, i) .and. a(cindex)<table%tbval(ccol, i+ddd) ) then
261 crow=i
262 cindex = cindex+1
263 dd = ddd
264 call gettablegrad( a, cindex, table, dd, crow, val1 )
265 crow = i+ddd
266 call gettablegrad( a, cindex, table, dd, crow, val2 )
267 lambda = (a(cindex-1)-table%tbval(ccol, i))/(table%tbval(ccol, crow)-table%tbval(ccol, i))
268 outa = (1.d0-lambda)*val1+ lambda* val2
269 exit
270 endif
271 enddo
272 endif
273 endif
274
275 end subroutine
276
279 subroutine fetch_tabledata( key, dict, outa, ierr, a )
280 character(len=*), intent(in) :: key
281 type(dict_struct), pointer :: dict
282 real(kind=kreal), intent(out) :: outa(:)
283 logical, intent(out) :: ierr
284 real(kind=kreal), intent(in), optional :: a(:)
285
286 type(dict_data), pointer :: dicval
287 integer :: nval, na, dd, crow, cindex
288
289 dicval => dict_get_key( dict, key )
290 ierr = .false.
291 if( .not. associated(dicval) ) then
292 ierr=.true.; return
293 endif
294 ! call print_table( dicval, 6 )
295
296 nval = dicval%tbcol-dicval%ndepends
297 if( nval /= size(outa) ) then
298 ierr=.true.; return
299 endif
300 if( dicval%tbrow==1 ) then
301 outa(:) =dicval%tbval(1:nval, 1); return
302 endif
303 if( .not. present(a) ) then
304 outa(:) =dicval%tbval(1:nval, 1); return
305 endif
306
307 cindex = 1; crow=1
308 dd = dicval%tbrow
309 na = 1
310 if( size(a) > dicval%ndepends ) then ! in case no dependent definition
311 na = size(a) - dicval%ndepends+1
312 endif
313 ! if( size(a) < dicval%ndepends ) then ! in case take no consider of dependent definition
314 ! na = dicval%ndepends- size(a)+1
315 ! cindex = na
316 ! endif
317 call gettabledata( a(na:), cindex, dicval, dd, crow, outa )
318 ! call finalize_table( dicval )
319
320 end subroutine
321
323 recursive subroutine gettabledata( a, cindex, table, dd, crow, outa )
324 real(kind=kreal), intent(in) :: a(:)
325 integer, intent(inout) :: cindex
326 type(dict_data) :: table
327 integer, intent(inout) :: dd, crow
328 real(kind=kreal), intent(out) :: outa(:)
329
330 integer :: i, ccol, ddd, nval
331 real(kind=kreal) :: lambda, val1(maxindex), val2(maxindex)
332
333 ddd = dd / table%tbindex(cindex)
334 ccol = table%tbcol-cindex+1
335 nval = size(outa)
336
337 if( ddd==1 ) then
338 if( a(cindex)<table%tbval(ccol, crow) ) then
339 outa(:) = table%tbval(1:ccol-1, crow)
340 elseif( a(cindex)>=table%tbval(ccol, crow+dd-1) ) then
341 outa(:) = table%tbval(1:ccol-1, crow+dd-1)
342 else
343 do i=crow, crow+dd-2
344 if( a(cindex)>=table%tbval(ccol, i) .and. a(cindex)<table%tbval(ccol, i+1) ) then
345 lambda = (a(cindex)-table%tbval(ccol, i))/(table%tbval(ccol, i+1)-table%tbval(ccol, i))
346 outa(:) = (1.d0-lambda)*table%tbval(1:ccol-1, i)+ lambda* table%tbval(1:ccol-1, i+1)
347 exit
348 endif
349 enddo
350 endif
351 else
352 if( a(cindex)<=table%tbval(ccol, crow) ) then
353 cindex = cindex+1
354 dd = ddd
355 call gettabledata( a, cindex, table, dd, crow, outa )
356 elseif( a(cindex)>=table%tbval(ccol, crow+dd-1) ) then
357 crow = crow+dd-ddd
358 cindex = cindex+1
359 dd = ddd
360 call gettabledata( a, cindex, table, dd, crow, outa )
361 else
362 do i=crow, crow+dd-2, ddd
363 if( a(cindex)==table%tbval(ccol, i) ) then
364 crow = i
365 cindex = cindex+1
366 dd = ddd
367 call gettabledata( a, cindex, table, dd, crow, outa )
368 exit
369 elseif( a(cindex)==table%tbval(ccol, i+ddd) ) then
370 crow = i+ddd
371 cindex = cindex+1
372 dd = ddd
373 call gettabledata( a, cindex, table, dd, crow, outa )
374 exit
375 elseif( a(cindex)>table%tbval(ccol, i) .and. a(cindex)<table%tbval(ccol, i+ddd) ) then
376 crow=i
377 cindex = cindex+1
378 dd = ddd
379 call gettabledata( a, cindex, table, dd, crow, val1(1:nval) )
380 crow = i+ddd
381 call gettabledata( a, cindex, table, dd, crow, val2(1:nval) )
382 lambda = (a(cindex-1)-table%tbval(ccol, i))/(table%tbval(ccol, crow)-table%tbval(ccol, i))
383 outa(:) = (1.d0-lambda)*val1(1:nval)+ lambda* val2(1:nval)
384 exit
385 endif
386 enddo
387 endif
388 endif
389
390 end subroutine
391
392
394 subroutine print_tabledata( dict, fname )
395 type(dict_struct), pointer :: dict
396 integer, intent(in) :: fname
397
398 type(linked_list), pointer :: current
399 integer :: i
400 do i = 1,size(dict%table)
401 if ( associated( dict%table(i)%list ) ) then
402 current => dict%table(i)%list
403 do while ( associated(current) )
404 if( trim(current%data%key) /= 'INIT' ) then
405 write( fname, * ) trim(current%data%key)
406 call print_table( current%data%value, fname )
407 endif
408 current => current%next
409 enddo
410 endif
411 enddo
412 end subroutine
413
414end module table_dicts
415
416
This module provides data structure table which would be dictionaried afterwards.
Definition: ttable.f90:7
subroutine finalize_table(table)
Definition: ttable.f90:81
integer, parameter maxnval
Definition: ttable.f90:16
subroutine print_table(table, fname)
Definition: ttable.f90:86
type(ttable), parameter dict_null
Definition: ttable.f90:30
subroutine tablecopy(lhs, rhs)
Definition: ttable.f90:115
logical function tablecompare(lhs, rhs)
Definition: ttable.f90:96
integer, parameter dict_key_length
Definition: ttable.f90:14
subroutine init_table(table, ndp, col, row, tbval)
Definition: ttable.f90:43
integer, parameter maxindex
Definition: ttable.f90:15
This module provides data structure of dictionaried table list.
Definition: ttable.f90:140
subroutine fetch_table(key, dict, dicval, ierr)
fetch a data table itself. P.A. it sould be deleted by users of this subroutine
Definition: ttable.f90:153
subroutine fetch_tabledata(key, dict, outa, ierr, a)
fetch a data by interpolation (This subroutine is specified for calculating temperature dependent har...
Definition: ttable.f90:280
subroutine fetch_tablegrad(key, a, dict, outa, ierr)
fetch a data by interpolation (This subroutine is specified for calculating temperature dependent har...
Definition: ttable.f90:183
subroutine print_tabledata(dict, fname)
Print our the contents of a dictionary.
Definition: ttable.f90:395
integer function fetch_tablerow(key, dict)
fetch a data table row
Definition: ttable.f90:169