1#define THIS_FILE "fdf.F90" 2!===================================================================== 3! 4! This file is part of the FDF package. 5! 6! This module implements an extended Fortran 90/95 interface 7! to the Flexible Data Format library of A. Garcia and J.M. Soler, 8! originally written in Fortran 77. 9! 10! FEATURES: 11! 12! a) Block pointers. 13! 14! Block content can be flexibly handled by means of a pointer 15! to a derived type 'block_fdf'. Typical usage: 16! 17! use fdf 18! type(block_fdf) :: bfdf 19! type(parsed_line), pointer :: pline 20! 21! if (fdf_block('SomeBlock', bfdf)) then 22! do while(fdf_bline(bfdf, pline)) 23! (process line 'integers|reals|values|names ...') 24! enddo 25! call fdf_bclose(bfdf) 26! endif 27! 28! The subroutine 'fdf_block' returns in 'bfdf' a structure used 29! to read the contents of the block. 30! 31! Routine fdf_bline returns in 'pline' the next non-blank parsed 32! line, non-comment line from the block, unless there are no more 33! lines, in which case it returns .FALSE. and 'pline' is undefined. 34! 35! Routine fdf_bclose runs the remaining lines in the block and ensures 36! the log may be used as input in subsequent entries. 37! 38! Routine 'backspace' moves the internal pointer of 'block_fdf' 39! structure to the previous line returned. 40! 41! Routine 'rewind' moves the internal pointer of 'block_fdf' structure 42! to the beginning of the block. 43! 44! b) Generic interface to scalar routines. 45! 46! The generic function 'fdf_get' can be used instead of any of the 47! scalar routines. The specific names are also accepted. 48! 49! c) Architecture support: this FDF implementation supports the following 50! architectures environments. 51! 52! 1) Thread-safe: The new implementation is thread-safe and will support 53! calling it from several OMP-threads executing in the same node. 54! 55! The implementation is as follows: fdf_init and fdf_shutdown are 56! SINGLE/CRITICAL sections that only one thread must execute. 57! On the other hand 'get'/'test' routines in FDF library are 58! thread-safe because each thread keeps its relative information 59! about the search/query that the caller program requests. 60! 61! 2) MPI-aware: For MPI executions, FDF library renames output/debugging 62! or log files to prevent overlaps of these files in a 63! shared/parallel filesystem. 64! 65! This option is enabled with _MPI_ macro, and superseded by 66! CLUSTER and BLOCKING macros. 67! 68! 3) Cluster filesystem: It is able to read the input FDF file from a 69! non-shared filesystem and broadcast the information to the rest 70! of the nodes in the execution. 71! 72! The implementation is as follows: the first node in the MPI rank 73! that is the owner of the input FDF file, process the file and 74! send it to the rest of the nodes in the MPI communicator. 75! 76! This option is enabled with CLUSTER macro. This option cannot be 77! used if BLOCKING macro is enabled. 78! 79! 4) Blocking reading: For huge executions (> 1.000 nodes) this option 80! is useful for shared/parallel filesystems where reading the same 81! file by several nodes could be a problem (collapsing system). 82! 83! The implementation is as follows: reading phase is done in blocking 84! pattern of size BLOCKSIZE nodes (configurable at compile time). 85! This means that the number of steps needed is STEPS = #NODES/BLOCKSIZE. 86! 87! This option is enabled with BLOCKING macro. This option cannot be 88! used if CLUSTER macro is enabled. 89! 90! Alberto Garcia, 1996-2007 91! Raul de la Cruz (BSC), September 2007 92! 93! 94!======================================================================== 95 96! AG: If running under MPI (flagged by the MPI symbol, choose the CLUSTER 97! mode of operation. 98! 99! 100#ifdef MPI 101# define _MPI_ 102# define CLUSTER 103# undef MPI 104#endif 105 106#ifndef _MPI_ 107# if defined(CLUSTER) || defined(BLOCKING) 108# define _MPI_ 109# endif 110#endif 111 112MODULE fdf 113 USE io_fdf 114 115 USE parse, only: parsed_line 116 USE parse, only: nintegers, nreals 117 USE parse, only: nvalues, nnames, ntokens 118 USE parse, only: integers, reals 119 USE parse, only: values, names, tokens, characters 120 USE parse, only: match 121 USE parse, only: digest, blocks, endblocks, labels 122 USE parse, only: destroy, setdebug, setlog, setmorphol 123 USE parse, only: nlists, lists 124 125 USE parse, only: search 126 USE parse, only: fdf_bsearch => search 127 USE parse, only: fdf_substring_search => substring_search 128 129 USE parse, only: serialize_pline, recreate_pline 130 USE parse, only: SERIALIZED_LENGTH 131 132 USE utils 133 USE prec 134 implicit none 135 136! User callable routines in FDF library 137 138! Start, stop and check parallelism in FDF system 139 public :: fdf_init, fdf_shutdown, fdf_parallel 140 141! Reading label functions 142 public :: fdf_get 143 public :: fdf_integer, fdf_single, fdf_double 144 public :: fdf_string, fdf_boolean 145 public :: fdf_physical, fdf_convfac 146 public :: fdf_islist, fdf_list 147 148! Returns the string associated with a mark line 149 public :: fdf_getline 150 151! Test if label is defined 152 public :: fdf_defined, fdf_isphysical, fdf_isblock 153 154! Allow to overwrite things in the FDF 155 public :: fdf_overwrite, fdf_removelabel, fdf_addline 156 157! Test if a label is used in obsolete or a deprecated state 158 public :: fdf_deprecated, fdf_obsolete 159 160! %block reading (processing each line) 161 public :: fdf_block, fdf_block_linecount 162 public :: fdf_bline, fdf_bbackspace, fdf_brewind, fdf_bclose 163 public :: fdf_bnintegers, fdf_bnreals, fdf_bnvalues, fdf_bnnames, fdf_bntokens 164 public :: fdf_bintegers, fdf_breals, fdf_bvalues, fdf_bnames, fdf_btokens 165 public :: fdf_bboolean, fdf_bphysical 166 public :: fdf_bnlists, fdf_blists 167 168! Match, search over blocks, and destroy block structure 169 public :: fdf_bmatch, fdf_bsearch, fdf_substring_search 170 171 public :: fdf_setoutput, fdf_setdebug 172 173! Private functions, non-callable 174 175! Main functions to build FDF structure (called in fdf_init) 176 private :: fdf_initdata, fdf_addtoken, fdf_readline 177 private :: fdf_read, fdf_readlabel, fdf_searchlabel 178 private :: fdf_open, fdf_close 179 180! Input/Output configuration 181 private :: fdf_input 182! private :: fdf_set_output_file 183 184! Destroy dynamic list of FDF structure (called in fdf_shutdown) 185 private :: fdf_destroy, fdf_destroy_dl 186 187! MPI init/finalize functions 188#ifdef _MPI_ 189 private :: fdf_mpi_init, fdf_mpi_finalize 190#endif 191 192! Reading functions for CLUSTER and BLOCKING configuration 193#ifdef CLUSTER 194 private :: setup_fdf_cluster, broadcast_fdf_struct 195#endif 196#ifdef BLOCKING 197 private :: fdf_readblocking 198#endif 199 200! Debugging functions, level and prints debugging info 201 private :: fdf_printfdf 202 203! Finds a label in the FDF herarchy 204 private :: fdf_locate 205 206! Dump function (for blocks) 207 private :: fdf_dump 208 209! Wrappers functions for block access, search, matching, 210! number and elements in the block (call to parse module) 211 interface fdf_bnintegers 212 module procedure nintegers 213 end interface 214 215 interface fdf_bnlists 216 module procedure nlists 217 end interface 218 219 interface fdf_bnreals 220 module procedure nreals 221 end interface 222 223 interface fdf_bnvalues 224 module procedure nvalues 225 end interface 226 227 interface fdf_bnnames 228 module procedure nnames 229 end interface 230 231 interface fdf_bntokens 232 module procedure ntokens 233 end interface 234 235 interface fdf_bintegers 236 module procedure integers 237 end interface 238 239 interface fdf_blists 240 module procedure lists 241 end interface 242 243 interface fdf_breals 244 module procedure reals 245 end interface 246 247 interface fdf_bvalues 248 module procedure values 249 end interface 250 251 interface fdf_bnames 252 module procedure names 253 end interface 254 255 interface fdf_btokens 256 module procedure tokens 257 end interface 258 259 interface fdf_bmatch 260 module procedure match 261 end interface 262 263! fdf_get wrapper for label functions 264 interface fdf_get 265 module procedure fdf_integer 266 module procedure fdf_single 267 module procedure fdf_double 268 module procedure fdf_boolean 269 module procedure fdf_string 270 module procedure fdf_physical 271 end interface 272 273 274! Unit numbers for input, output, error notification, and 275! debugging output (the latter active if fdf_debug is true) 276 logical, private :: fdf_debug = .FALSE., & 277 fdf_debug2 = .FALSE., & 278 fdf_started = .FALSE., & 279 fdf_output = .FALSE. 280 281 integer(ip), parameter, private :: maxdepth = 7 282 integer(ip), parameter, private :: maxFileNameLength = 300 283 integer(ip), private :: ndepth 284 integer(ip), private :: fdf_in(maxdepth) 285 integer(ip), private :: fdf_out, fdf_err, fdf_log 286 287! MPI variables (id, number of MPI tasks) 288#ifdef _MPI_ 289 logical, private :: mpiflag 290 integer(ip), private :: rank, ntasks 291#endif 292 293! Structure for searching inside fdf blocks 294 type, public :: block_fdf 295 character(len=MAX_LENGTH) :: label 296 type(line_dlist), pointer :: mark => null() 297 end type block_fdf 298 299! Dynamic list for parsed_line structures 300 type, public :: line_dlist 301 character(len=MAX_LENGTH) :: str 302 type(parsed_line), pointer :: pline => null() 303 type(line_dlist), pointer :: next => null() 304 type(line_dlist), pointer :: prev => null() 305 end type line_dlist 306 307! FDF data structure (first and last lines) 308 type, private :: fdf_file 309 integer(ip) :: nlines 310 type(line_dlist), pointer :: first => null() 311 type(line_dlist), pointer :: last => null() 312 end type fdf_file 313 314! Input FDF file 315 type(fdf_file), pointer, private :: file_in => null() 316 317 318 319! Define by default all the others inherit module entities as privated 320! avoiding redefinitions of entities in several module files with same name 321 public :: parsed_line ! Structure for searching inside fdf blocks 322 public :: leqi ! For legacy support (old codes) 323 private 324 325 326 CONTAINS 327 328! 329! Initialization for fdf. 330! 331 SUBROUTINE fdf_init( fileInput, fileOutput, unitInput ) 332 implicit none 333!------------------------------------------------------------- Input Variables 334 character(len=*),optional,intent(in):: fileInput, fileOutput 335 integer, optional,intent(in):: unitInput 336 337#ifndef FDF_DEBUG 338!------------------------------------------------------------- Local Variables 339 integer(ip) :: debug_level, output_level 340#endif 341 character(len=256) :: filedebug 342 character(len=maxFileNameLength):: filein, fileout 343 344!----------------------------------------------------------------------- BEGIN 345!$OMP SINGLE 346 ! Prevent the user from opening two head files 347 if (fdf_started) then 348 call die('FDF module: fdf_init', 'Head file already set', & 349 THIS_FILE, __LINE__, fdf_err) 350 endif 351 352#ifdef _MPI_ 353 call fdf_mpi_init() 354#endif 355 call fdf_initdata() 356 357 call io_geterr(fdf_err) 358 359 ! Set in/out file names, if fileInput and fileOutput are not present 360 call set_file_names( filein, fileout, & 361 fileInput, fileOutput, unitInput ) 362 filedebug = trim(fileout) // ".debug" 363 364#ifdef FDF_DEBUG 365 ! To monitor the parsing and the build-up of the 366 ! fdf data structures in all nodes 367 call fdf_setdebug(2,filedebug) 368 call fdf_setoutput(2,fileout) ! All nodes print output 369#endif 370 371!! call fdf_set_output_file(fileout) 372 373 call fdf_input(filein) 374 375 fdf_started = .TRUE. 376 377#ifndef FDF_DEBUG 378 ! Flags within the fdf file itself. 379 380 ! At this point only the final fdf data structure will be shown, 381 ! for level >= 2 382 debug_level = fdf_get('fdf-debug', 0) 383 call fdf_setdebug(debug_level,filedebug) 384 385 ! The default is to have output only in the master node 386 output_level = fdf_get('fdf-output', 1) 387 call fdf_setoutput(output_level,fileout) 388#endif 389 390 if (fdf_debug2) call fdf_printfdf() 391 392#ifdef _MPI_ 393 call fdf_mpi_finalize() 394#endif 395!$OMP END SINGLE 396 RETURN 397!------------------------------------------------------------------------- END 398 END SUBROUTINE fdf_init 399 400 401 SUBROUTINE set_file_names( fileIn, fileOut, & 402 optFileIn, optFileOut, unitIn ) 403 ! If present, copies input arguments optFileIn/Out to fileIn/Out. 404 ! If absent, generates In/Out file names. If unitIn is present, and it is 405 ! a named file, returns it as fileIn. If not, it copies input to a new 406 ! file and returns its name. If .not.present(unitIn) => unitIn=5. 407 ! If optFileIn is present, unitIn is ignored. 408 implicit none 409 character(len=*),intent(out):: & 410 fileIn, &! Name of file to be used as input 411 fileOut ! Name of file to be used as output 412 character(len=*),optional,intent(in):: & 413 optFileIn, &! Optional argument with input file name 414 optFileOut ! Optional argument with output file name 415 integer,optional,intent(in):: & 416 unitIn ! Optional input file unit (not used if present(optFileIn)) 417 418 integer:: count, ierr, iostat, iu, iuIn 419 logical:: named, opened 420 character(len=MAX_LENGTH*2) line 421 character(len=maxFileNameLength) fileName 422 423!------------------------------------------------------------------------- BEGIN 424#ifdef _MPI_ 425 if (rank==0) then 426#endif 427 428 ! Find a job-specific number 429 call system_clock( count ) 430 count = mod(count,100000) 431 432 ! Set output file name 433 if (present(optFileOut)) then 434 if (len(trim(optFileOut)) > len(fileOut)) & 435 call die('FDF module: set_file_names', & 436 'Parameter maxFileNameLength too small.' // & 437 'Terminating.', THIS_FILE, __LINE__, fdf_err, rc=ierr) 438 fileOut = optFileOut 439 else ! set a job-specific file name 440 write(fileOut,'(a,i5.5,a)') 'fdf_',count,'.log' 441 endif 442 443 ! Set input file 444 if (present(optFileIn)) then ! just copy the file name 445 if (len(trim(optFileIn)) > len(fileIn)) & 446 call die('FDF module: set_file_names', & 447 'Parameter maxFileNameLength too small.' // & 448 'Terminating.', THIS_FILE, __LINE__, fdf_err, rc=ierr) 449 fileIn = optFileIn 450 else ! find or set a file name 451 452 ! Find input file unit 453 if (present(unitIn)) then ! use given unit (possibly 5) 454 iuIn = unitIn 455 else ! assume standard input 456 iuIn = 5 457 endif 458 459 ! Find file name associated with given unit 460 if (iuIn==5) then ! no valid file name 461 fileName = ' ' 462 else ! check if this is a named file 463 inquire(unit=iuIn,opened=opened) 464 if (opened) then 465 inquire(unit=iuIn,named=named) 466 if (named) then ! inquire file name 467 inquire(unit=iuIn,name=fileName) 468 else ! no valid file name 469 fileName = ' ' 470 endif ! (named) 471 else 472 call die('FDF module: set_file_names', 'Input unit not opened.' // & 473 'Terminating.', THIS_FILE, __LINE__, fdf_err, rc=ierr) 474 endif ! (opened) 475 endif ! (iuIn==5) 476 477 ! Set input file name, possibly after copying input to it 478 if (fileName==' ') then ! not a valid file 479 write(fileIn,'(a,i5.5,a)') & 480 'INPUT_TMP_',count,'.fdf' ! new file's name 481 call io_assign(iu) ! new file's unit 482 open(iu,file=trim(fileIn),form='formatted') ! open new file 483 do 484 read(iuIn,iostat=iostat,fmt='(a)') line ! read line from old unit 485 if (iostat/=0 ) exit 486 write(iu,'(a)') trim(line) ! write line to new file 487 enddo 488 call io_close(iu) ! close new file 489 else ! valid file 490 fileIn = fileName 491 endif ! (fileName=='stdin') 492 493 endif ! (present(optFileIn)) 494 495#ifdef _MPI_ 496 endif ! (rank==0) 497#endif 498!--------------------------------------------------------------------------- END 499 END SUBROUTINE set_file_names 500 501! 502! Initialize MPI subsystem if the application calling/using FDF 503! library is not running with MPI enabled. 504! 505#ifdef _MPI_ 506 SUBROUTINE fdf_mpi_init() 507 use mpi_siesta 508 implicit none 509!--------------------------------------------------------------- Local Variables 510 integer(ip) :: ierr 511 512!------------------------------------------------------------------------- BEGIN 513 call MPI_Initialized(mpiflag, ierr) 514 if (ierr .ne. MPI_SUCCESS) then 515 call die('FDF module: fdf_mpi_init', 'Error initializing MPI system.' // & 516 'Terminating.', THIS_FILE, __LINE__, fdf_err, rc=ierr) 517 endif 518 519 if (.not. mpiflag) then 520 call MPI_Init(ierr) 521 if (ierr .ne. MPI_SUCCESS) then 522 call die('FDF module: fdf_mpi_init', 'Error initializing MPI system.' // & 523 'Terminating.', THIS_FILE, __LINE__, fdf_err, rc=ierr) 524 endif 525 endif 526 527 call MPI_Comm_Rank(MPI_COMM_WORLD, rank, ierr) 528 if (ierr .ne. MPI_SUCCESS) then 529 call die('FDF module: fdf_mpi_init', 'Error getting MPI comm rank.' // & 530 'Terminating.', THIS_FILE, __LINE__, fdf_err, rc= ierr) 531 endif 532 533 call MPI_Comm_Size(MPI_COMM_WORLD, ntasks, ierr) 534 if (ierr .ne. MPI_SUCCESS) then 535 call die('FDF module: fdf_mpi_init', 'Error getting MPI comm size.' // & 536 'Terminating.', THIS_FILE, __LINE__, fdf_err, rc=ierr) 537 endif 538 RETURN 539!--------------------------------------------------------------------------- END 540 END SUBROUTINE fdf_mpi_init 541#endif 542 543! 544! Finalize MPI subsystem if the application calling/using FDF 545! library is not running with MPI enabled. 546! 547#ifdef _MPI_ 548 SUBROUTINE fdf_mpi_finalize() 549 use mpi_siesta 550 implicit none 551!--------------------------------------------------------------- Local Variables 552 integer(ip) :: ierr 553 554!------------------------------------------------------------------------- BEGIN 555 if (.not. mpiflag) then 556 call MPI_Finalize(ierr) 557 if (ierr .ne. MPI_SUCCESS) then 558 call die('FDF module: fdf_mpi_finalize', 'Error finalizing MPI system.' // & 559 'Terminating.', THIS_FILE, __LINE__, fdf_err, rc=ierr) 560 endif 561 endif 562 563 RETURN 564!--------------------------------------------------------------------------- END 565 END SUBROUTINE fdf_mpi_finalize 566#endif 567 568! 569! Reads the input file depending on the configuration of the system: 570! Shared Filesystem, Cluster Filesystem or Blocking input access 571! 572 SUBROUTINE fdf_input(filein) 573 implicit none 574!--------------------------------------------------------------- Input Variables 575 character(*) :: filein 576 577!------------------------------------------------------------------------- BEGIN 578#ifdef CLUSTER 579!! call fdf_readcluster(filein) 580 call setup_fdf_cluster(filein) 581#elif defined(BLOCKING) 582 call fdf_readblocking(filein) 583#else 584 call fdf_read(filein) 585#endif 586 587 if (fdf_output) write(fdf_out,'(a,a,a,i3)') '#FDF module: Opened ', filein, & 588 ' for input. Unit:', fdf_in(1) 589 590 RETURN 591!--------------------------------------------------------------------------- END 592 END SUBROUTINE fdf_input 593 594! 595! Reading code for Cluster Architecture, where the filesystem 596! is not shared along all the nodes in the system. 597! 598#ifdef CLUSTER 599 SUBROUTINE fdf_readcluster(filein) 600 use mpi_siesta 601 implicit none 602!--------------------------------------------------------------- Input Variables 603 character(*) :: filein 604 605!--------------------------------------------------------------- Local Variables 606 character(80) :: msg 607 character(256) :: fileinTmp 608 integer(ip) :: ierr, texist_send, texist_recv 609#ifdef SOPHISTICATED_SEARCH 610 logical :: file_exist 611#endif 612!------------------------------------------------------------------------- BEGIN 613! Tests if the running node has the input file: 614! If found: texist_send = rank 615! Else : texist_send = error_code (ntasks + 1) 616 617#ifdef SOPHISTICATED_SEARCH 618 INQUIRE(file=filein, exist=file_exist) 619 if (file_exist) then 620 texist_send = rank 621 else 622 texist_send = ntasks + 1 623 endif 624 625 call MPI_AllReduce(texist_send, texist_recv, 1, MPI_INTEGER, & 626 MPI_MIN, MPI_COMM_WORLD, ierr) 627 if (ierr .ne. MPI_SUCCESS) then 628 call die('FDF module: fdf_readcluster', 'Error in MPI_AllReduce (task_exist).' // & 629 'Terminating.', THIS_FILE, __LINE__, fdf_err, rc=ierr) 630 endif 631#else 632! 633! Simplify: Assume node 0 has the file 634! 635 texist_recv = 0 636#endif 637 638! The node owner of the input file send the data to the other ones 639! if none node has an input file with such name abort the application 640 if (texist_recv .eq. ntasks + 1) then 641 write(msg,*) 'No node found in the cluster with ', & 642 ' input file ', filein,'. Terminating.' 643 call die('FDF module: fdf_readcluster', msg, & 644 THIS_FILE, __LINE__, fdf_err, rc=ierr) 645 else 646 if (texist_recv .eq. rank) then 647 call fdf_read(filein) 648 call fdf_sendInput() 649!!Debug call MPI_Barrier( MPI_COMM_WORLD, ierr ) 650 if (fdf_output) write(fdf_out,*) '#FDF module: Node', rank, 'reading/sending', & 651 ' input file ', filein 652 else 653 call fdf_recvInput(texist_recv, filein, fileinTmp) 654!!Debug call MPI_Barrier( MPI_COMM_WORLD, ierr ) 655 call fdf_read(fileinTmp) 656 if (fdf_output) write(fdf_out,*) '#FDF module: Node', rank, 'receiving input', & 657 ' file from', texist_recv, 'to ', TRIM(fileinTmp) 658 endif 659 endif 660 RETURN 661!--------------------------------------------------------------------------- END 662 END SUBROUTINE fdf_readcluster 663! 664! 665 SUBROUTINE setup_fdf_cluster(filein) 666! 667! A more efficient alternative to fdf_sendInput/fdf_recvInput 668! that avoids the creation of scratch files by non-root nodes. 669! 670! Alberto Garcia, April 2011 671 672 use mpi_siesta 673 implicit none 674!--------------------------------------------------------------- Input Variables 675 character(*), intent(in) :: filein ! File name 676 677!--------------------------------------------------------------- Local Variables 678 character(80) :: msg 679 character(256) :: fileinTmp 680 integer(ip) :: ierr, texist_send, reading_node 681#ifdef SOPHISTICATED_SEARCH 682 logical :: file_exist 683#endif 684!------------------------------------------------------------------------- BEGIN 685! Tests if the running node has the input file: 686! If found: texist_send = rank 687! Else : texist_send = error_code (ntasks + 1) 688 689#ifdef SOPHISTICATED_SEARCH 690 INQUIRE(file=filein, exist=file_exist) 691 if (file_exist) then 692 texist_send = rank 693 else 694 texist_send = ntasks + 1 695 endif 696 697 call MPI_AllReduce(texist_send, reading_node, 1, MPI_INTEGER, & 698 MPI_MIN, MPI_COMM_WORLD, ierr) 699 if (ierr .ne. MPI_SUCCESS) then 700 call die('FDF module: fdf_readcluster', 'Error in MPI_AllReduce (task_exist).' // & 701 'Terminating.', THIS_FILE, __LINE__, fdf_err, rc=ierr) 702 endif 703#else 704! 705! Simplify: Assume node 0 has the file 706! 707 reading_node = 0 708#endif 709 710! If no node has an input file with such name abort the application 711 if (reading_node .eq. ntasks + 1) then 712 write(msg,*) 'No node found in the cluster with ', & 713 ' input file ', filein,'. Terminating.' 714 call die('FDF module: fdf_readcluster', msg, & 715 THIS_FILE, __LINE__, fdf_err, rc=ierr) 716 endif 717 718! The root node reads and digests the input file 719! and sends the serialized fdf_file data structure to the other nodes 720 721 if (rank == reading_node) then 722 call fdf_read(filein) 723 endif 724 call broadcast_fdf_struct(reading_node) 725 726 RETURN 727!--------------------------------------------------------------------------- END 728 END SUBROUTINE setup_fdf_cluster 729 730! 731! Broadcast complete fdf structure 732! 733 SUBROUTINE broadcast_fdf_struct(reading_node) 734 use mpi_siesta 735 implicit none 736 737 integer, intent(in) :: reading_node ! Node which contains the struct 738 739!--------------------------------------------------------------- Local Variables 740 character, pointer :: bufferFDF(:) => null() 741 integer(ip) :: i, j, k, ierr, nlines 742 743!------------------------------------------------------------------------- BEGIN 744 745 if (rank == reading_node) then 746 nlines = file_in%nlines 747 endif 748 749 call MPI_Bcast(nlines, 1, & 750 MPI_INTEGER, reading_node, MPI_COMM_WORLD, ierr) 751 if (ierr .ne. MPI_SUCCESS) then 752 call die('FDF module: broadcast_fdf', 'Error Broadcasting nlines.' // & 753 'Terminating.', THIS_FILE, __LINE__, fdf_err, rc=ierr) 754 endif 755 756 ALLOCATE(bufferFDF(nlines*SERIALIZED_LENGTH), stat=ierr) 757 if (ierr .ne. 0) then 758 call die('FDF module: broadcast_fdf', 'Error allocating bufferFDF', & 759 THIS_FILE, __LINE__, fdf_err, rc=ierr) 760 endif 761 762 if (rank == reading_node) then 763 call serialize_fdf_struct(bufferFDF) 764 endif 765 766 call MPI_Bcast(bufferFDF, size(bufferFDF), & 767 MPI_CHARACTER, reading_node, MPI_COMM_WORLD, ierr) 768 if (ierr .ne. MPI_SUCCESS) then 769 call die('FDF module: broadcast_fdf', 'Error Broadcasting bufferFDF.' // & 770 'Terminating.', THIS_FILE, __LINE__, fdf_err, rc=ierr) 771 endif 772 if (rank /= reading_node) then 773 call recreate_fdf_struct(nlines,bufferFDF) 774 endif 775 776 DEALLOCATE(bufferFDF) 777 778 RETURN 779!--------------------------------------------------------------------------- END 780 END SUBROUTINE broadcast_fdf_struct 781#endif 782! 783! Reading code for Blocking read. The reading of the input file is 784! splitted in several steps of size BLOCKSIZE. The number of steps 785! needed to read the file in all the nodes depends on the total nodes, 786! being ntasks/BLOCKINGSIZE. 787! 788! With this method we can avoid a colapse of a global parallel filesystem 789! if the execution of the application runs over a huge amount of nodes. 790! 791#ifdef BLOCKING 792# ifndef BLOCKSIZE 793# define BLOCKSIZE 2 794# endif 795 SUBROUTINE fdf_readblocking(filein) 796 use mpi_siesta 797 implicit none 798!--------------------------------------------------------------- Input Variables 799 character(*) :: filein 800 801!--------------------------------------------------------------- Local Variables 802 integer(ip) :: i, ierr 803 804!------------------------------------------------------------------------- BEGIN 805 806 do i= 0, ntasks-1, BLOCKSIZE 807 if ((rank .ge. i) .and. (rank .le. i+BLOCKSIZE-1)) then 808 call fdf_read(filein) 809 if (fdf_output) write(fdf_out,*) '#FDF module: Task', rank, 'reading input', & 810 ' file ', filein, ' in step', (i/BLOCKSIZE)+1 811 endif 812 813 call MPI_Barrier(MPI_COMM_WORLD, ierr) 814 if (ierr .ne. MPI_SUCCESS) then 815 call die('FDF module: fdf_readblocking', 'Error in MPI_Barrier (fdf_read).' // & 816 'Terminating.', THIS_FILE, __LINE__, fdf_err, rc=ierr) 817 endif 818 enddo 819 820 RETURN 821!--------------------------------------------------------------------------- END 822 END SUBROUTINE fdf_readblocking 823#endif 824 825! 826! Read an input file (and include files) and builds memory 827! structure that will contain the data and will help in searching 828! 829 RECURSIVE SUBROUTINE fdf_read(filein, blocklabel) 830 implicit none 831!--------------------------------------------------------------- Input Variables 832 character(*) :: filein 833 character(*), optional :: blocklabel 834 835!--------------------------------------------------------------- Local Variables 836 logical :: dump 837 logical, allocatable :: found(:) 838 character(80) :: msg 839 character(len=MAX_LENGTH) :: label, inc_file 840 character(len=MAX_LENGTH*2):: line 841 integer(ip) :: i, ierr, ntok, ind_less, nlstart 842 type(parsed_line), pointer :: pline 843 844!------------------------------------------------------------------------- BEGIN 845! Open reading input file 846 call fdf_open(filein) 847 848! Read each input data line 849 if (PRESENT(blocklabel)) then 850 label = blocklabel 851 else 852 label = ' ' 853 endif 854 do while (fdf_readline(line)) 855 856! Check if valid data (tokens, non-blank) 857 pline => digest(line) 858 ntok = ntokens(pline) 859 if (ntok .ne. 0) then 860 861! Find different special cases in the input files 862! (%block, %endblock, %include, Label1 Label2 ... < Filename) 863 864! %block directive 865 ind_less = search('<', pline) 866 if (search('%block', pline) .eq. 1) then 867 868! No label found in %block directive 869 if (ntok .eq. 1) then 870 write(msg,*) '%block label not found in ', TRIM(filein) 871 call die('FDF module: fdf_read', msg, & 872 THIS_FILE, __LINE__, fdf_err) 873 endif 874 875! %block Label < Filename [ %dump ] 876 if (ind_less .eq. 3) then 877 878 if (ntok .ge. 4) then 879! Test if %dump is present 880 if (search('%dump', pline) .eq. 5) then 881 dump = .TRUE. 882 else 883 dump = .FALSE. 884 endif 885 886! Add begin, body and end sections of block 887 label = tokens(pline, 2) 888 inc_file = tokens(pline, 4) 889 call destroy(pline) 890 line = '%block ' // label 891 pline => digest(line) 892 call setmorphol(1, 'b', pline) 893 call setmorphol(2, 'l', pline) 894 call fdf_addtoken(line, pline) 895 nullify(pline) ! it is stored in line 896 897 nlstart = file_in%nlines 898 call fdf_read(inc_file, label) 899 900! Warn if block 'label' is empty 901 if ((nlstart - file_in%nlines) .eq. 0) then 902 write(msg,*) 'FDF module: fdf_read: block ', & 903 TRIM(label), ' is empty...' 904 call warn(msg) 905 endif 906 907 line = '%endblock ' // label 908 pline => digest(line) 909 call setmorphol(1, 'e', pline) 910 call setmorphol(2, 'l', pline) 911 call fdf_addtoken(line, pline) 912 nullify(pline) ! it is stored in line 913 914! Dump included file to fileout 915 if (dump) call fdf_dump(label) 916 label = ' ' 917 918! Filename not found in %block directive 919 else 920 write(msg,*) '%block filename not found in ', TRIM(filein) 921 call die('FDF module: fdf_read', msg, & 922 THIS_FILE, __LINE__, fdf_err) 923 endif 924 925! %block Label 926 elseif (ind_less .eq. -1) then 927 label = tokens(pline, 2) 928 call setmorphol(1, 'b', pline) 929 call setmorphol(2, 'l', pline) 930 call fdf_addtoken(line, pline) 931 nullify(pline) ! it is stored in line 932 nlstart = file_in%nlines 933 934! Bad format in %block directive 935 else 936 write(msg,*) 'Bad ''<'' %block format in ', TRIM(filein) 937 call die('FDF module: fdf_read', msg, & 938 THIS_FILE, __LINE__, fdf_err) 939 endif 940 941! %endblock directive 942 elseif (search('%endblock', pline) .eq. 1) then 943! Check if %block exists before %endblock 944 if (label .eq. ' ') then 945 write(msg,*) 'Bad %endblock found in ', TRIM(filein) 946 call die('FDF module: fdf_read', msg, & 947 THIS_FILE, __LINE__, fdf_err) 948 else 949! Warn if block 'label' is empty 950 if ((nlstart - file_in%nlines) .eq. 0) then 951 write(msg,*) 'FDF module: fdf_read: block ', & 952 TRIM(label), ' is empty...' 953 call warn(msg) 954 endif 955 956 call destroy(pline) 957 line = '%endblock ' // label 958 pline => digest(line) 959 call setmorphol(1, 'e', pline) 960 call setmorphol(2, 'l', pline) 961 call fdf_addtoken(line, pline) 962 nullify(pline) ! it is stored in line 963 label = ' ' 964 endif 965 966! %include Filename directive 967 elseif (search('%include', pline) .eq. 1) then 968! Check if include filename is specified 969 if (ntok .eq. 1) then 970 write(msg,*) 'Filename on %include not found in ', TRIM(filein) 971 call die('FDF module: fdf_read', msg, & 972 THIS_FILE, __LINE__, fdf_err) 973 else 974 inc_file = tokens(pline, 2) 975 call fdf_read(inc_file) 976 endif 977 978 ! Clean pline (we simply insert the next file) 979 call destroy(pline) 980 981! Label1 Label2 ... < Filename directive 982 elseif (ind_less .ne. -1) then 983! Check if '<' is in a valid position 984 if (ind_less .eq. 1) then 985 write(msg,*) 'Bad ''<'' found in ', TRIM(filein) 986 call die('FDF module: fdf_read', msg, & 987 THIS_FILE, __LINE__, fdf_err) 988 989! Check if '<' filename is specified 990 elseif (ind_less .eq. ntok) then 991 write(msg,*) 'Filename not found after ''<'' in ', TRIM(filein) 992 call die('FDF module: fdf_read', msg, & 993 THIS_FILE, __LINE__, fdf_err) 994 else 995! Search label(s) in Filename 996 inc_file = tokens(pline, ind_less+1) 997 ALLOCATE(found(ind_less-1), stat=ierr) 998 if (ierr .ne. 0) then 999 call die('FDF module: fdf_read', 'Error allocating found', & 1000 THIS_FILE, __LINE__, fdf_err, rc=ierr) 1001 endif 1002 1003! If label(s) not found in such Filename throw an error 1004 found = .FALSE. 1005 if (.not. fdf_readlabel(ind_less-1, pline, & 1006 inc_file, found)) then 1007 i = 1 1008 do while ((i .le. ind_less-1) .and. (found(i))) 1009 i = i + 1 1010 enddo 1011 label = tokens(pline, i) 1012 write(msg,*) 'Label ', TRIM(label), & 1013 ' not found in ', TRIM(inc_file) 1014 call die('FDF module: fdf_read', msg, & 1015 THIS_FILE, __LINE__, fdf_err) 1016 endif 1017 1018 call destroy(pline) 1019 DEALLOCATE(found) 1020 endif 1021 1022! Add remaining kind of tokens to dynamic list as labels 1023 else 1024 if (label .eq. ' ') call setmorphol(1, 'l', pline) 1025 call fdf_addtoken(line, pline) 1026 nullify(pline) ! it is stored in line 1027 endif 1028 else 1029! Destroy parsed_line structure if no elements 1030 call destroy(pline) 1031 endif 1032 enddo 1033 1034! Close one level of input file 1035 if ((.not. PRESENT(blocklabel)) .and. (label .ne. ' ')) then 1036 write(msg,*) '%endblock ', TRIM(label), & 1037 ' not found in ', TRIM(filein) 1038 call die('FDF module: fdf_read', msg, THIS_FILE, __LINE__, fdf_err) 1039 endif 1040 call fdf_close() 1041 1042 RETURN 1043!--------------------------------------------------------------------------- END 1044 END SUBROUTINE fdf_read 1045 1046! 1047! Read an input file (and include files) searching labels to 1048! include them in memory structure that will contain the data 1049! 1050 RECURSIVE FUNCTION fdf_readlabel(nelem, plabel, filein, found) result(readlabel) 1051 implicit none 1052!--------------------------------------------------------------- Input Variables 1053 character(*) :: filein 1054 integer(ip) :: nelem 1055 type(parsed_line), pointer :: plabel 1056 1057!-------------------------------------------------------------- Output Variables 1058 logical :: readlabel 1059 logical :: found(nelem) 1060 1061!--------------------------------------------------------------- Local Variables 1062 logical :: dump, found_elem 1063 logical, pointer :: found_loc(:) 1064 character(80) :: msg 1065 character(len=MAX_LENGTH*2):: line 1066 character(len=MAX_LENGTH) :: label, inc_file 1067 integer(ip) :: i, ierr, ntok, ind_less, nlstart 1068 integer(ip) :: elem, nelem_loc 1069 integer(ip), pointer :: found_index(:) 1070 type(parsed_line), pointer :: pline 1071 1072!------------------------------------------------------------------------- BEGIN 1073! Open input file with labels 1074 call fdf_open(filein) 1075 1076! While not reach to end of file and found all labels 1077 do while (fdf_readline(line) .and. (.not. ALL(found))) 1078 1079! Check if valid data (tokens, non-blank) 1080 pline => digest(line) 1081 ntok = ntokens(pline) 1082 if (ntok .ne. 0) then 1083 1084! Find different special cases in the input files 1085! (%block, %endblock, %include, Label1 Label2 ... < Filename) 1086 1087! %block directive 1088 ind_less = search('<', pline) 1089 if (search('%block', pline) .eq. 1) then 1090 1091! No label found in %block directive 1092 if (ntok .eq. 1) then 1093 write(msg,*) '%block label not found in ', TRIM(filein) 1094 call die('FDF module: fdf_readlabel', msg, & 1095 THIS_FILE, __LINE__, fdf_err) 1096 endif 1097 1098! %block Label < Filename [ %dump ] 1099 if (ind_less .eq. 3) then 1100 1101 if (ntok .ge. 4) then 1102! Test if %dump is present 1103 if (search('%dump', pline) .eq. 5) then 1104 dump = .TRUE. 1105 else 1106 dump = .FALSE. 1107 endif 1108 1109 label = tokens(pline, 2) 1110 elem = fdf_searchlabel(found, nelem, label, plabel) 1111 1112 inc_file = tokens(pline, 4) 1113 call destroy(pline) 1114 1115! If match with any label add [begin, body, end] of block 1116 if (elem .ne. -1) then 1117 line = '%block ' // label 1118 pline => digest(line) 1119 call setmorphol(1, 'b', pline) 1120 call setmorphol(2, 'l', pline) 1121 call fdf_addtoken(line, pline) 1122 1123 nlstart = file_in%nlines 1124 call fdf_read(inc_file, label) 1125 1126! Warn if block 'label' is empty 1127 if ((nlstart - file_in%nlines) .eq. 0) then 1128 write(msg,*) 'FDF module: fdf_readlabel: block ', & 1129 TRIM(label), ' is empty...' 1130 call warn(msg) 1131 endif 1132 1133 line = '%endblock ' // label 1134 pline => digest(line) 1135 call setmorphol(1, 'e', pline) 1136 call setmorphol(2, 'l', pline) 1137 call fdf_addtoken(line, pline) 1138 1139! Dump included file to fileout 1140 if (dump) call fdf_dump(label) 1141 1142 found(elem) = .TRUE. 1143 label = ' ' 1144 endif 1145 1146! Filename not found in %block directive 1147 else 1148 write(msg,*) 'Filename on %block not found in ', TRIM(filein) 1149 call die('FDF module: fdf_readlabel', msg, & 1150 THIS_FILE, __LINE__, fdf_err) 1151 endif 1152 1153! %block Label 1154 elseif (ind_less .eq. -1) then 1155 label = tokens(pline, 2) 1156 elem = fdf_searchlabel(found, nelem, label, plabel) 1157 found_elem = .TRUE. 1158 1159! If match with any label add [begin,body,end] of block 1160 if (elem .ne. -1) then 1161 call setmorphol(1, 'b', pline) 1162 call setmorphol(2, 'l', pline) 1163 call fdf_addtoken(line, pline) 1164 nlstart = file_in%nlines 1165 1166 found_elem = .FALSE. 1167 do while (fdf_readline(line) .and. (.not. found_elem)) 1168 pline => digest(line) 1169 if (ntokens(pline) .ne. 0) then 1170 if (search('%endblock', pline) .eq. 1) then 1171! Warn if block 'label' is empty 1172 if ((nlstart - file_in%nlines) .eq. 0) then 1173 write(msg,*) 'FDF module: fdf_readlabel: block ', & 1174 TRIM(label), ' is empty...' 1175 call warn(msg) 1176 endif 1177 1178 call destroy(pline) 1179 line = '%endblock ' // label 1180 pline => digest(line) 1181 call setmorphol(1, 'e', pline) 1182 call setmorphol(2, 'l', pline) 1183 label = ' ' 1184 1185 found_elem = .TRUE. 1186 found(elem) = .TRUE. 1187 endif 1188 call fdf_addtoken(line, pline) 1189 endif 1190 enddo 1191 1192! Move to the end of the block 1193 else 1194 call destroy(pline) 1195 1196 found_elem = .FALSE. 1197 do while (fdf_readline(line) .and. (.not. found_elem)) 1198 pline => digest(line) 1199 if (search('%endblock', pline) .eq. 1) then 1200 label = ' ' 1201 found_elem = .TRUE. 1202 endif 1203 call destroy(pline) 1204 enddo 1205 endif 1206 1207! Error due to %endblock not found 1208 if (.not. found_elem) then 1209 write(msg,*) '%endblock ', TRIM(label), & 1210 ' not found in ', TRIM(filein) 1211 call die('FDF module: fdf_readlabel', msg, & 1212 THIS_FILE, __LINE__, fdf_err) 1213 endif 1214 1215! Bad format in %block directive 1216 else 1217 write(msg,*) 'Bad ''<'' %block format in ', TRIM(filein) 1218 call die('FDF module: fdf_readlabel', msg, & 1219 THIS_FILE, __LINE__, fdf_err) 1220 endif 1221 1222! %endblock directive 1223 elseif (search('%endblock', pline) .eq. 1) then 1224! Bad if %endblock exists before %block 1225 write(msg,*) 'Bad %endblock found in ', TRIM(filein) 1226 call die('FDF module: fdf_readlabel', msg, & 1227 THIS_FILE, __LINE__, fdf_err) 1228 1229! %include Filename directive 1230 elseif (search('%include', pline) .eq. 1) then 1231! Check if include filename is specified 1232 if (ntok .eq. 1) then 1233 write(msg,*) 'Filename on %include not found in ', TRIM(filein) 1234 call die('FDF module: fdf_readlabel', msg, & 1235 THIS_FILE, __LINE__, fdf_err) 1236 else 1237 inc_file = tokens(pline, 2) 1238 call destroy(pline) 1239 readlabel = fdf_readlabel(nelem, plabel, inc_file, found) 1240 endif 1241 1242! Label1 Label2 ... < Filename directive 1243 elseif (ind_less .ne. -1) then 1244! Check if '<' is in a valid position 1245 if (ind_less .eq. 1) then 1246 write(msg,*) 'Bad ''<'' found in ', TRIM(filein) 1247 call die('FDF module: fdf_readlabel', msg, & 1248 THIS_FILE, __LINE__, fdf_err) 1249 1250! Check if '<' filename is specified 1251 elseif (ind_less .eq. ntok) then 1252 write(msg,*) 'Filename not found after ''<'' in ', TRIM(filein) 1253 call die('FDF module: fdf_readlabel', msg, & 1254 THIS_FILE, __LINE__, fdf_err) 1255 else 1256! Search label(s) in Filename 1257 line = ' ' 1258 nelem_loc = 0 1259 ALLOCATE(found_index(ind_less-1), stat=ierr) 1260 if (ierr .ne. 0) then 1261 call die('FDF module: fdf_readlabel', 'Error allocating found_index', & 1262 THIS_FILE, __LINE__, fdf_err, rc=ierr) 1263 endif 1264 do i= 1, ind_less-1 1265 label = tokens(pline, i) 1266 elem = fdf_searchlabel(found, nelem, label, plabel) 1267 if (elem .ne. -1) then 1268 line = TRIM(line) // ' ' // TRIM(label) 1269 nelem_loc = nelem_loc + 1 1270 found_index(nelem_loc) = elem 1271 endif 1272 enddo 1273 1274! Process Filename if any label found 1275 if (nelem_loc .ge. 1) then 1276 inc_file = tokens(pline, ind_less+1) 1277 call destroy(pline) 1278 1279 ALLOCATE(found_loc(nelem_loc), stat=ierr) 1280 if (ierr .ne. 0) then 1281 call die('FDF module: fdf_readlabel', 'Error allocating found_loc', & 1282 THIS_FILE, __LINE__, fdf_err, rc=ierr) 1283 endif 1284 1285 found_loc = .FALSE. 1286 1287! If label(s) not found in such Filename throw an error 1288 pline => digest(line) 1289 if (.not. fdf_readlabel(nelem_loc, pline, & 1290 inc_file, found_loc)) then 1291 i = 1 1292 do while ((i .le. nelem_loc) .and. (found_loc(i))) 1293 i = i + 1 1294 enddo 1295 label = tokens(pline, i) 1296 write(msg,*) 'Label ', TRIM(label), ' not found in ', TRIM(inc_file) 1297 call die('FDF module: fdf_readlabel', msg, & 1298 THIS_FILE, __LINE__, fdf_err) 1299 else 1300! Merge results if all labels found 1301 do i= 1, nelem_loc 1302 found(found_index(i)) = found_loc(i) 1303 enddo 1304 endif 1305 1306 DEALLOCATE(found_index) 1307 endif 1308 1309 DEALLOCATE(found_loc) 1310 call destroy(pline) 1311 endif 1312 1313! Label [ Value ] directive 1314 else 1315 elem = fdf_searchlabel(found, nelem, tokens(pline, 1), plabel) 1316 1317! If match with any label add it 1318 if (elem .ne. -1) then 1319 call setmorphol(1, 'l', pline) 1320 call fdf_addtoken(line, pline) 1321 found(elem) = .TRUE. 1322 else 1323! Destroy parsed_line structure if no label found 1324 call destroy(pline) 1325 endif 1326 endif 1327 1328 else 1329! Destroy parsed_line structure if no label found 1330 call destroy(pline) 1331 endif 1332 enddo 1333 1334! Close input file with labels 1335 call fdf_close() 1336 1337 readlabel = ALL(found) 1338 RETURN 1339!--------------------------------------------------------------------------- END 1340 END FUNCTION fdf_readlabel 1341 1342! 1343! Search a label in a set of not 'found' tokens given by plabel and nelem. 1344! Returns the index on plabel of the token that matches with label. 1345! 1346 FUNCTION fdf_searchlabel(found, nelem, label, plabel) 1347 implicit none 1348!--------------------------------------------------------------- Input Variables 1349 integer(ip) :: nelem 1350 logical :: found(nelem) 1351 character(*) :: label 1352 type(parsed_line), pointer :: plabel 1353 1354!-------------------------------------------------------------- Output Variables 1355 integer(ip) :: fdf_searchlabel 1356 1357!--------------------------------------------------------------- Local Variables 1358 logical :: found_elem 1359 integer(ip) :: i 1360 1361!------------------------------------------------------------------------- BEGIN 1362 i = 1 1363 found_elem = .FALSE. 1364 fdf_searchlabel = -1 1365 do while ((i .le. nelem) .and. (.not. found_elem)) 1366 1367 if (.not. found(i)) then 1368 if (labeleq(label, tokens(plabel, i))) then 1369 found_elem = .TRUE. 1370 fdf_searchlabel = i 1371 endif 1372 endif 1373 i = i + 1 1374 enddo 1375 1376 RETURN 1377!--------------------------------------------------------------------------- END 1378 END FUNCTION fdf_searchlabel 1379 1380! 1381! Dumps the content of a block 1382! 1383 SUBROUTINE fdf_dump(label) 1384 implicit none 1385!--------------------------------------------------------------- Input Variables 1386 character(*) :: label 1387 1388!--------------------------------------------------------------- Local Variables 1389 character(80) :: msg 1390 type(block_fdf) :: bfdf 1391 type(parsed_line), pointer :: pline 1392 1393!------------------------------------------------------------------------- BEGIN 1394 fdf_started = .TRUE. 1395 1396 if (.not. fdf_block(label, bfdf)) then 1397 write(msg,*) 'block ', label, 'to dump not found' 1398 call die('FDF module: fdf_dump', msg, THIS_FILE, __LINE__, fdf_err) 1399 endif 1400 1401! fdf_bline prints each block line in fdf_out 1402 do while(fdf_bline(bfdf, pline)) 1403 enddo 1404 1405 fdf_started = .FALSE. 1406 1407 RETURN 1408!--------------------------------------------------------------------------- END 1409 END SUBROUTINE fdf_dump 1410 1411! 1412! Init FDF file structure 1413! 1414 SUBROUTINE fdf_initdata() 1415 implicit none 1416!--------------------------------------------------------------- Local Variables 1417 integer(ip) :: ierr 1418 1419!------------------------------------------------------------------------- BEGIN 1420 ndepth = 0 1421 1422 ALLOCATE(file_in, stat=ierr) 1423 if (ierr .ne. 0) then 1424 call die('FDF module: fdf_initdata', 'Error allocating file_in', & 1425 THIS_FILE, __LINE__, fdf_err, rc=ierr) 1426 endif 1427 1428 file_in%nlines = 0 1429 NULLIFY(file_in%first) 1430 NULLIFY(file_in%last) 1431 1432 RETURN 1433!--------------------------------------------------------------------------- END 1434 END SUBROUTINE fdf_initdata 1435 1436! 1437! Add a line individually to the dynamic list of parsed lines 1438! This can not include block's and is restricted to key values 1439! 1440 SUBROUTINE fdf_addline(line) 1441 implicit none 1442!--------------------------------------------------------------- Input Variables 1443 character(len=*) :: line 1444 1445!--------------------------------------------------------------- Local Variables 1446 integer(ip) :: ntok 1447 type(parsed_line), pointer :: pline 1448 1449!------------------------------------------------------------------------- BEGIN 1450 1451! Check if valid data (tokens, non-blank) 1452 pline => digest(line) 1453 1454 call setmorphol(1, 'l', pline) 1455 call fdf_addtoken(line, pline) 1456 1457 if (fdf_debug2) then 1458 write(fdf_log,*) '***FDF_ADDLINE********************************' 1459 write(fdf_log,*) 'Line:', TRIM(line) 1460 write(fdf_log,*) '**********************************************' 1461 endif 1462 1463 END SUBROUTINE fdf_addline 1464 1465! 1466! Remove a line from the dynamic list of parsed lines 1467! This can not include block's and is restricted to key values 1468! 1469 SUBROUTINE fdf_removelabel(label) 1470 implicit none 1471!--------------------------------------------------------------- Input Variables 1472 character(len=*) :: label 1473 1474!--------------------------------------------------------------- Local Variables 1475 type(line_dlist), pointer :: mark 1476 1477!------------------------------------------------------------------------- BEGIN 1478 1479 do while ( fdf_locate(label,mark) ) 1480 1481 if (fdf_debug2) then 1482 write(fdf_log,*) '***FDF_REMOVELABEL*******************************' 1483 write(fdf_log,*) 'Line:', TRIM(mark%str) 1484 write(fdf_log,*) 'Label:', trim(label) 1485 write(fdf_log,*) '**********************************************' 1486 endif 1487 1488 ! To circumvent the first/last line in the fdf-file 1489 ! we have to check for the existance of the 1490 ! first/last mark being the one removed. 1491 ! That special case *must* correct the first/last 1492 ! tokens. 1493 if ( associated(mark,target=file_in%first) ) then 1494 file_in%first => mark%next 1495 end if 1496 if ( associated(mark,target=file_in%last) ) then 1497 file_in%last => mark%prev 1498 end if 1499 1500 ! Remove the label from the dynamic list 1501 call destroy(mark%pline) 1502 if ( associated(mark%prev) ) then 1503 mark%prev%next => mark%next 1504 end if 1505 if ( associated(mark%next) ) then 1506 mark%next%prev => mark%prev 1507 end if 1508 DEALLOCATE(mark) 1509 1510 NULLIFY(mark) 1511 end do 1512 1513 END SUBROUTINE fdf_removelabel 1514 1515! 1516! Overwrite label line in dynamic list of parsed lines 1517! 1518 SUBROUTINE fdf_overwrite(line) 1519!--------------------------------------------------------------- Input Variables 1520 character(len=*) :: line 1521 1522!--------------------------------------------------------------- Local Variables 1523 type(parsed_line), pointer :: pline 1524 character(len=MAX_LENGTH) :: label 1525 1526 integer :: ierr 1527 1528 pline => digest(line) 1529 if ( search('%block', pline) == 1 .or. & 1530 search('%endblock', pline) == 1 ) then 1531 1532 ! We do not allow this in a single line 1533 call die('FDF module: fdf_overwrite', 'Error overwriting block (not implemented)', & 1534 THIS_FILE, __LINE__, fdf_err, rc=ierr) 1535 1536 else if ( search('%include', pline) == 1 ) then 1537 1538 ! We do not allow this in a single line 1539 call die('FDF module: fdf_overwrite', 'Error overwriting flags from input file (not implemented)', & 1540 THIS_FILE, __LINE__, fdf_err, rc=ierr) 1541 1542 else if ( search('<', pline) /= -1 ) then 1543 1544 ! We do not allow this in a single line 1545 call die('FDF module: fdf_overwrite', 'Error piping in overwriting (not implemented)', & 1546 THIS_FILE, __LINE__, fdf_err, rc=ierr) 1547 1548 else 1549 1550 label = tokens(pline,1) 1551 call setmorphol(1, 'l', pline) 1552 call fdf_removelabel(label) 1553 1554 ! Add token to the list of fdf-flags 1555 ! Since we add it directly we shouldn't destroy the pline 1556 call fdf_addtoken(line, pline) 1557 if ( fdf_debug ) then 1558 write(fdf_log,'(2a)') '---> Overwriting token: ', trim(label) 1559 end if 1560 1561 end if 1562 1563 END SUBROUTINE fdf_overwrite 1564 1565! 1566! Add a token to the dynamic list of parsed lines 1567! 1568 SUBROUTINE fdf_addtoken(line, pline) 1569 implicit none 1570!--------------------------------------------------------------- Input Variables 1571 character(len=*) :: line 1572 type(parsed_line), pointer :: pline 1573 1574!--------------------------------------------------------------- Local Variables 1575 integer(ip) :: i, ierr 1576 type(line_dlist), pointer :: mark 1577 1578!------------------------------------------------------------------------- BEGIN 1579 ALLOCATE(mark, stat=ierr) 1580 if (ierr .ne. 0) then 1581 call die('FDF module: fdf_addtoken', 'Error allocating mark', & 1582 THIS_FILE, __LINE__, fdf_err, rc=ierr) 1583 endif 1584 1585 mark%str = line 1586 mark%pline => pline 1587 NULLIFY(mark%next) 1588 1589 if (ASSOCIATED(file_in%first)) then 1590 mark%prev => file_in%last 1591 file_in%last%next => mark 1592 else 1593 NULLIFY(mark%prev) 1594 file_in%first => mark 1595 endif 1596 1597 file_in%last => mark 1598 file_in%nlines = file_in%nlines + 1 1599 1600 if (fdf_debug2) then 1601 write(fdf_log,*) '***FDF_ADDTOKEN*******************************' 1602 write(fdf_log,*) 'Line:', TRIM(mark%str) 1603 write(fdf_log,*) 'Ntokens:', mark%pline%ntokens 1604 do i= 1, mark%pline%ntokens 1605 write(fdf_log,*) ' Token:', trim(tokens(pline,i)), & 1606 ' (', mark%pline%id(i), ')' 1607 enddo 1608 write(fdf_log,*) '**********************************************' 1609 endif 1610 1611 RETURN 1612!--------------------------------------------------------------------------- END 1613 END SUBROUTINE fdf_addtoken 1614 1615! 1616! Opens a file for FDF processing. 1617! 1618 SUBROUTINE fdf_open(filename) 1619 implicit none 1620!-------------------------------------------------------------- Output Variables 1621 character(*) :: filename 1622 1623!--------------------------------------------------------------- Local Variables 1624 logical :: file_exists 1625 character(80) :: msg 1626 integer(ip) :: lun 1627 1628!------------------------------------------------------------------------- BEGIN 1629 ndepth = ndepth + 1 1630 if (ndepth .gt. maxdepth) then 1631 call die('FDF module: fdf_open', 'Too many nested fdf files...', & 1632 THIS_FILE, __LINE__, fdf_err) 1633 endif 1634 1635 if (leqi(filename, 'stdin')) then 1636 lun = INPUT_UNIT 1637 if (fdf_debug) write(fdf_log,'(a,i1,a)') & 1638 '---> Reading from standard input [DEPTH:', ndepth,'] ' 1639 else 1640 call io_assign(lun) 1641 1642 INQUIRE(file=filename, exist=file_exists) 1643 if (file_exists) then 1644 open(unit=lun, file=filename, status='old', form='formatted') 1645 REWIND(lun) 1646 if (fdf_debug) write(fdf_log,'(a,i1,a,a)') & 1647 '---> Opened [DEPTH:', ndepth,'] ', TRIM(filename) 1648 else 1649 write(msg,'(a,a)') 'Cannot open ', TRIM(filename) 1650 call die('FDF module: fdf_open', msg, THIS_FILE, __LINE__, fdf_err) 1651 endif 1652 endif 1653 1654 fdf_in(ndepth) = lun 1655 REWIND(fdf_in(ndepth)) 1656 1657 RETURN 1658!--------------------------------------------------------------------------- END 1659 END SUBROUTINE fdf_open 1660 1661! 1662! Closes currently opened fdf file 1663! 1664 SUBROUTINE fdf_close() 1665 implicit none 1666!------------------------------------------------------------------------- BEGIN 1667 if (ndepth .ge. 1) then 1668 call io_close(fdf_in(ndepth)) 1669 if (fdf_debug) & 1670 write(fdf_log,'(a,i1,a)') '---> Closed [DEPTH:', ndepth,']' 1671 ndepth = ndepth - 1 1672 endif 1673 1674 RETURN 1675!--------------------------------------------------------------------------- END 1676 END SUBROUTINE fdf_close 1677 1678! 1679! Set output file for FDF subsystem 1680! 1681 SUBROUTINE fdf_set_output_file(fileout) 1682 implicit none 1683!----------------------------------------------------- Input Variables 1684 character(len=*), intent(in) :: fileout 1685 1686!----------------------------------------------------- Local Variables 1687 character(256) :: fileouttmp 1688!----------------------------------------------------- BEGIN 1689 call io_assign(fdf_out) 1690 1691#ifdef _MPI_ 1692 if (rank /= 0) then 1693 fileouttmp = "/tmp/" // trim(fileout) // "." // i2s(rank) 1694 else 1695 fileouttmp = fileout 1696 endif 1697#else 1698 fileouttmp = fileout 1699#endif 1700 1701#ifdef FDF_DEBUG 1702 ! 1703 ! If debugging, all the nodes use named log files 1704 ! 1705 open( unit=fdf_out, file=TRIM(fileouttmp), form='formatted', & 1706 access='sequential', status='replace' ) 1707#else 1708 ! 1709 ! Only the master node opens a named log file in the current dir. 1710 ! Non-master nodes use a scratch file. 1711 ! These log files tend to be quite small, so there 1712 ! should not be problems such as filling up filesystems. 1713 ! ... your mileage might vary. This is a grey area. 1714 ! Some compilers allow the user to specify where scratch 1715 ! files go. 1716 ! 1717#ifdef _MPI_ 1718 if (rank /= 0) then 1719 open( unit=fdf_out, form='formatted', & 1720 access='sequential', status='scratch' ) 1721 else 1722 open( unit=fdf_out, file=TRIM(fileouttmp), form='formatted', & 1723 access='sequential', status='replace' ) 1724 end if 1725#else 1726 open( unit=fdf_out, file=TRIM(fileouttmp), form='formatted', & 1727 access='sequential', status='replace' ) 1728#endif 1729#endif 1730 1731 RETURN 1732!--------------------------------------------------------------------------- END 1733 END SUBROUTINE fdf_set_output_file 1734 1735! 1736! Frees and shutdown FDF system 1737! 1738 SUBROUTINE fdf_shutdown() 1739 implicit none 1740!------------------------------------------------------------------------- BEGIN 1741!$OMP SINGLE 1742 if (fdf_started) then 1743 call fdf_destroy(file_in) 1744 fdf_started = .FALSE. 1745 1746 call io_close(fdf_out) 1747 endif 1748!$OMP END SINGLE 1749 1750 RETURN 1751!--------------------------------------------------------------------------- END 1752 END SUBROUTINE fdf_shutdown 1753 1754! 1755! Destroy the fdf_file structure for the input file 1756! 1757 SUBROUTINE fdf_destroy(fdfp) 1758 implicit none 1759!-------------------------------------------------------------- Output Variables 1760 type(fdf_file), pointer :: fdfp 1761 1762!------------------------------------------------------------------------- BEGIN 1763 if (ASSOCIATED(fdfp%first)) call fdf_destroy_dl(fdfp%first) 1764 DEALLOCATE(fdfp) 1765 1766 RETURN 1767!--------------------------------------------------------------------------- END 1768 END SUBROUTINE fdf_destroy 1769 1770! 1771! Destroy recursively the dynamic list of parsed lines 1772! 1773 RECURSIVE SUBROUTINE fdf_destroy_dl(dlp) 1774 implicit none 1775!-------------------------------------------------------------- Output Variables 1776 type(line_dlist), pointer :: dlp 1777 1778!------------------------------------------------------------------------- BEGIN 1779 if (ASSOCIATED(dlp%next)) call fdf_destroy_dl(dlp%next) 1780 call destroy(dlp%pline) 1781 DEALLOCATE(dlp) 1782 1783 RETURN 1784!--------------------------------------------------------------------------- END 1785 END SUBROUTINE fdf_destroy_dl 1786 1787! 1788! Tells when FDF library has been compiled with parallel 1789! execution support. Depends on _MPI_ macro. 1790! 1791 FUNCTION fdf_parallel() 1792 implicit none 1793!-------------------------------------------------------------- Output Variables 1794 logical :: fdf_parallel 1795 1796!------------------------------------------------------------------------- BEGIN 1797#ifdef _MPI_ 1798 fdf_parallel = .TRUE. 1799#else 1800 fdf_parallel = .FALSE. 1801#endif 1802 1803 RETURN 1804!--------------------------------------------------------------------------- END 1805 END FUNCTION fdf_parallel 1806 1807! 1808! Read a line of the 'ndepth' input file, returning .TRUE. if 1809! there are more lines to read from input file, .FALSE. otherwise. 1810! 1811 FUNCTION fdf_readline(line) 1812 implicit none 1813!-------------------------------------------------------------- Output Variables 1814 logical :: fdf_readline 1815 character(*) :: line 1816 1817!--------------------------------------------------------------- Local Variables 1818 integer(ip) :: stat 1819 1820!------------------------------------------------------------------------- BEGIN 1821 read(fdf_in(ndepth), '(a)', iostat=stat) line 1822 1823 if (stat .eq. 0) then 1824 fdf_readline = .TRUE. 1825 if (fdf_debug2) write(fdf_log, '(a,a76)') 'fdf_readline > ', line 1826 else 1827 fdf_readline = .FALSE. 1828 endif 1829 1830 RETURN 1831!--------------------------------------------------------------------------- END 1832 END FUNCTION fdf_readline 1833 1834! 1835! Returns the string line associated with the mark pointer 1836! in the FDF herarchy. If mark is not associated returns ''. 1837! 1838 FUNCTION fdf_getline(mark) 1839 implicit none 1840!--------------------------------------------------------------- Input Variables 1841 type(line_dlist), pointer :: mark 1842 1843!-------------------------------------------------------------- Output Variables 1844 character(len=MAX_LENGTH) :: fdf_getline 1845 1846!------------------------------------------------------------------------- BEGIN 1847 if (ASSOCIATED(mark)) then 1848 fdf_getline = mark%str 1849 else 1850 fdf_getline = ' ' 1851 endif 1852 1853 RETURN 1854!--------------------------------------------------------------------------- END 1855 END FUNCTION fdf_getline 1856 1857! 1858! Prints all the fdf_file structure of the input file(s) 1859! 1860 SUBROUTINE fdf_printfdf() 1861 implicit none 1862!--------------------------------------------------------------- Local Variables 1863 integer(ip) :: i, ntokens 1864 character*1 :: id 1865 type(line_dlist), pointer :: dlp 1866 character(len=MAX_LENGTH) :: tok 1867 1868!------------------------------------------------------------------------- BEGIN 1869 dlp => file_in%first 1870 1871 write(fdf_log,*) '*** FDF Memory Structure Summary: ************' 1872 do while (ASSOCIATED(dlp)) 1873 ntokens = dlp%pline%ntokens 1874 write(fdf_log,*) 'Line:', TRIM(dlp%str) 1875 write(fdf_log,*) 'Ntokens:', ntokens 1876 do i= 1, ntokens 1877 tok = tokens(dlp%pline,i) 1878 id = dlp%pline%id(i) 1879 write(fdf_log,*) ' Token:', trim(tok), '(', dlp%pline%id(i), ')' 1880 enddo 1881 dlp => dlp%next 1882 enddo 1883 write(fdf_log,*) '**********************************************' 1884 1885 RETURN 1886!--------------------------------------------------------------------------- END 1887 END SUBROUTINE fdf_printfdf 1888 1889! 1890! Returns an integer associated with label 'label', or the default 1891! value if label is not found in the fdf file. 1892! Optionally can return a pointer to the line found. 1893! 1894 FUNCTION fdf_integer(label, default, line) 1895 implicit none 1896!--------------------------------------------------------------- Input Variables 1897 character(*) :: label 1898 integer(ip) :: default 1899 1900!-------------------------------------------------------------- Output Variables 1901 integer(ip) :: fdf_integer 1902 type(line_dlist), pointer, optional :: line 1903 1904!--------------------------------------------------------------- Local Variables 1905 character(80) :: msg 1906 type(line_dlist), pointer :: mark 1907 1908!------------------------------------------------------------------------- BEGIN 1909! Prevents using FDF routines without initialize 1910 if (.not. fdf_started) then 1911 call die('FDF module: fdf_integer', 'FDF subsystem not initialized', & 1912 THIS_FILE, __LINE__, fdf_err) 1913 endif 1914 1915 if (fdf_locate(label, mark)) then 1916 if (.not. match(mark%pline, 'li')) then 1917 write(msg,*) 'no integer value for ', label 1918 call die('FDF module: fdf_integer', msg, THIS_FILE, __LINE__, fdf_err) 1919 endif 1920 1921 fdf_integer = integers(mark%pline, 1, 1) 1922 if (fdf_output) write(fdf_out,'(a,5x,i10)') label, fdf_integer 1923 else 1924 fdf_integer = default 1925 if (fdf_output) write(fdf_out,'(a,i10,5x,a)') label, default, '# default value' 1926 endif 1927 1928 if (PRESENT(line)) line = mark 1929 1930 RETURN 1931!--------------------------------------------------------------------------- END 1932 END FUNCTION fdf_integer 1933 1934! 1935! Returns true or false whether or not the label 'label' is 1936! a value with units or not. 1937! I.e. it returns true if the line has the form lvn, if not found, or not lvn, 1938! it returns false. 1939! 1940 FUNCTION fdf_isphysical(label) 1941 implicit none 1942!--------------------------------------------------------------- Input Variables 1943 character(*) :: label 1944 1945!-------------------------------------------------------------- Output Variables 1946 logical :: fdf_isphysical 1947 1948!--------------------------------------------------------------- Local Variables 1949 type(line_dlist), pointer :: mark 1950 1951!------------------------------------------------------------------------- BEGIN 1952! Prevents using FDF routines without initialize 1953 if (.not. fdf_started) then 1954 call die('FDF module: fdf_isphysical', 'FDF subsystem not initialized', & 1955 THIS_FILE, __LINE__, fdf_err) 1956 endif 1957 1958 if (fdf_locate(label, mark)) then 1959 fdf_isphysical = match(mark%pline, 'lvn') 1960 else 1961 fdf_isphysical = .false. 1962 endif 1963 if (fdf_output) write(fdf_out,'(a,5x,l10)') "#:physical? " // label, fdf_isphysical 1964 1965 RETURN 1966!--------------------------------------------------------------------------- END 1967 END FUNCTION fdf_isphysical 1968 1969! 1970! Returns true or false whether or not the label 'label' is 1971! a list or not, you cannot get the line out from this routine 1972! 1973 FUNCTION fdf_islist(label) 1974 implicit none 1975!--------------------------------------------------------------- Input Variables 1976 character(*) :: label 1977 1978!-------------------------------------------------------------- Output Variables 1979 logical :: fdf_islist 1980 1981!--------------------------------------------------------------- Local Variables 1982 type(line_dlist), pointer :: mark 1983 1984!------------------------------------------------------------------------- BEGIN 1985! Prevents using FDF routines without initialize 1986 if (.not. fdf_started) then 1987 call die('FDF module: fdf_islist', 'FDF subsystem not initialized', & 1988 THIS_FILE, __LINE__, fdf_err) 1989 endif 1990 1991 if (fdf_locate(label, mark)) then 1992 ! if it is a list: 1993 fdf_islist = match(mark%pline, 'la') 1994 else 1995 fdf_islist = .false. 1996 endif 1997 if (fdf_output) write(fdf_out,'(a,5x,l10)') "#:list? " // label, fdf_islist 1998 1999 RETURN 2000!--------------------------------------------------------------------------- END 2001 END FUNCTION fdf_islist 2002 2003! 2004! Returns a list with label 'label', or the default 2005! value if label is not found in the fdf file. 2006! 2007 SUBROUTINE fdf_list(label,ni,list,line) 2008 implicit none 2009!--------------------------------------------------------------- Input Variables 2010 character(*) :: label 2011 integer(ip) :: ni 2012 2013!-------------------------------------------------------------- Output Variables 2014 integer(ip) :: list(ni) 2015 type(line_dlist), pointer, optional :: line 2016 2017!--------------------------------------------------------------- Local Variables 2018 character(80) :: msg 2019 type(line_dlist), pointer :: mark 2020 integer(ip) :: lni, llist(1) 2021 2022!------------------------------------------------------------------------- BEGIN 2023! Prevents using FDF routines without initialize 2024 if (.not. fdf_started) then 2025 call die('FDF module: fdf_list', 'FDF subsystem not initialized', & 2026 THIS_FILE, __LINE__, fdf_err) 2027 endif 2028 2029 if (fdf_locate(label, mark)) then 2030 if (.not. match(mark%pline, 'la')) then 2031 write(msg,*) 'no list value for ', label 2032 call die('FDF module: fdf_list', msg, THIS_FILE, __LINE__, fdf_err) 2033 endif 2034 2035 ! Retrieve length of list 2036 lni = -1 2037 call lists(mark%pline,1,lni,llist) 2038 if ( ni <= 0 ) then 2039 ! the user has requested size... 2040 ni = lni 2041 else 2042 ! the list is not long enough 2043 if ( ni < lni ) then 2044 write(msg, '(2a,2(a,i0))')'List ', trim(label), & 2045 ' container too small: ', ni, ' versus ', lni 2046 call die('FDF module: fdf_list', trim(msg), & 2047 THIS_FILE, __LINE__, fdf_err) 2048 end if 2049 call lists(mark%pline,1,ni,list) 2050 end if 2051 2052 ! find a way to write out the list anyway 2053 !if (fdf_output) write(fdf_out,'(a,5x,i10)') label, lni 2054 else 2055 write(msg,*) 'no list value for ', label 2056 call die('FDF module: fdf_list', msg, THIS_FILE, __LINE__, fdf_err) 2057 endif 2058 2059 if (PRESENT(line)) line = mark 2060 2061 RETURN 2062!--------------------------------------------------------------------------- END 2063 END SUBROUTINE fdf_list 2064 2065! 2066! Returns a string associated with label 'label', or the default 2067! string if label is not found in the fdf file. 2068! Optionally can return a pointer to the line found. 2069! 2070 FUNCTION fdf_string(label, default, line) 2071 implicit none 2072!--------------------------------------------------------------- Input Variables 2073 character(*) :: label 2074 character(*) :: default 2075 2076!-------------------------------------------------------------- Output Variables 2077 character(80) :: fdf_string 2078 type(line_dlist), pointer, optional :: line 2079 2080!--------------------------------------------------------------- Local Variables 2081 type(line_dlist), pointer :: mark 2082 2083!------------------------------------------------------------------------- BEGIN 2084! Prevents using FDF routines without initialize 2085 if (.not. fdf_started) then 2086 call die('FDF module: fdf_string', 'FDF subsystem not initialized', & 2087 THIS_FILE, __LINE__, fdf_err) 2088 endif 2089 2090 if (fdf_locate(label, mark)) then 2091 if (ntokens(mark%pline) < 2) then 2092 fdf_string = "" 2093 if (fdf_output) write(fdf_out,'(a,5x,a)') label, & 2094 "# *** Set to empty string *** " 2095 else 2096 ! Get all the characters spanning the space from the second to 2097 ! the last token 2098 fdf_string = characters(mark%pline, ind_init=2, ind_final=-1) 2099 if (fdf_output) write(fdf_out,'(a,5x,a)') label, fdf_string 2100 endif 2101 else 2102 fdf_string = default 2103 if (fdf_output) write(fdf_out,'(a,5x,a,5x,a)') label, default, '# default value' 2104 endif 2105 2106 if (PRESENT(line)) line = mark 2107 2108 RETURN 2109!--------------------------------------------------------------------------- END 2110 END FUNCTION fdf_string 2111 2112! 2113! Returns true if label 'label' appears by itself or in the form 2114! label {yes,true,.true.,t,y} (case insensitive). 2115! 2116! Returns false if label 'label' appears in the form 2117! label {no,false,.false.,f,n} (case insensitive). 2118! 2119! If label is not found in the fdf file, fdf_boolean returns the 2120! LOGICAL variable default. 2121! 2122! Optionally can return a pointer to the line found. 2123! 2124 FUNCTION fdf_boolean(label, default, line) 2125 implicit none 2126!--------------------------------------------------------------- Input Variables 2127 character(*) :: label 2128 logical :: default 2129 2130!-------------------------------------------------------------- Output Variables 2131 logical :: fdf_boolean 2132 type(line_dlist), pointer, optional :: line 2133 2134!--------------------------------------------------------------- Local Variables 2135 character(80) :: msg, valstr 2136 type(line_dlist), pointer :: mark 2137 2138 2139!------------------------------------------------------------------------- BEGIN 2140! Prevents using FDF routines without initialize 2141 if (.not. fdf_started) then 2142 call die('FDF module: fdf_boolean', 'FDF subsystem not initialized', & 2143 THIS_FILE, __LINE__, fdf_err) 2144 endif 2145 2146 if (fdf_locate(label, mark)) then 2147 2148! If the label appears by itself, we interpret it as .true. 2149 if (ntokens(mark%pline) .ne. 1) then 2150 2151! Look for second word 2152 valstr = names(mark%pline, 1, 1) 2153 2154 if (is_true(valstr)) then 2155 fdf_boolean = .TRUE. 2156 if (fdf_output) write(fdf_out,'(a,5x,l10)') label, fdf_boolean 2157 2158 elseif (is_false(valstr)) then 2159 fdf_boolean = .FALSE. 2160 if (fdf_output) write(fdf_out,'(a,5x,l10)') label, fdf_boolean 2161 2162 else 2163 write(msg,*) 'unexpected logical value ', label, ' = ', valstr 2164 call die('FDF module: fdf_boolean', msg, & 2165 THIS_FILE, __LINE__, fdf_err) 2166 endif 2167 else 2168 fdf_boolean = .TRUE. 2169 if (fdf_output) write(fdf_out,'(a,5x,l10,5x,a)') label, fdf_boolean, & 2170 '# label by itself' 2171 endif 2172 else 2173 fdf_boolean = default 2174 if (fdf_output) write(fdf_out,'(a,5x,l10,5x,a)') label, default, '# default value' 2175 endif 2176 2177 if (PRESENT(line)) line = mark 2178 2179 RETURN 2180 2181 CONTAINS 2182 2183 logical function is_true(valstr) result(a) 2184 character(len=*), intent(in) :: valstr 2185 a = leqi(valstr, 'yes') .or. leqi(valstr, 'true') .or. & 2186 leqi(valstr, '.true.') .or. leqi(valstr, 't') .or. & 2187 leqi(valstr, 'y') 2188 end function is_true 2189 2190 logical function is_false(valstr) result(a) 2191 character(len=*), intent(in) :: valstr 2192 a = leqi(valstr, 'no') .or. leqi(valstr, 'false') .or. & 2193 leqi(valstr, '.false.') .or. leqi(valstr, 'f') .or. & 2194 leqi(valstr, 'n') 2195 end function is_false 2196 2197!--------------------------------------------------------------------------- END 2198 END FUNCTION fdf_boolean 2199 2200! 2201! Returns true if label 'label' appears by itself or in the form 2202! label {yes,true,.true.,t,y} (case insensitive). 2203! 2204! Returns false if label 'label' appears in the form 2205! label {no,false,.false.,f,n} (case insensitive). 2206! 2207! If label is not found in the fdf file, fdf_boolean returns the 2208! LOGICAL variable default. 2209! 2210! Optionally can return a pointer to the line found. 2211! 2212 FUNCTION fdf_bboolean(pline,ind,after) 2213 implicit none 2214!--------------------------------------------------------------- Input Variables 2215 integer(ip), intent(in) :: ind 2216 integer(ip), intent(in), optional :: after 2217 type(parsed_line), pointer :: pline 2218 2219!-------------------------------------------------------------- Output Variables 2220 logical :: fdf_bboolean 2221 2222!--------------------------------------------------------------- Local Variables 2223 character(80) :: msg, valstr 2224 type(line_dlist), pointer :: mark 2225 2226 2227!------------------------------------------------------------------------- BEGIN 2228! Prevents using FDF routines without initialize 2229 if (.not. fdf_started) then 2230 call die('FDF module: fdf_bboolean', 'FDF subsystem not initialized', & 2231 THIS_FILE, __LINE__, fdf_err) 2232 endif 2233 2234 if (ind <= nnames(pline,after=after)) then 2235 2236 valstr = names(pline,ind,after=after) 2237 2238 if (is_true(valstr)) then 2239 fdf_bboolean = .TRUE. 2240 if (fdf_output) write(fdf_out,'(a,5x,l10)') valstr, fdf_bboolean 2241 2242 elseif (is_false(valstr)) then 2243 fdf_bboolean = .FALSE. 2244 if (fdf_output) write(fdf_out,'(a,5x,l10)') valstr, fdf_bboolean 2245 2246 else 2247 write(msg,*) 'unexpected logical value ', valstr 2248 call die('FDF module: fdf_bboolean', msg, & 2249 THIS_FILE, __LINE__, fdf_err) 2250 endif 2251 else 2252 fdf_bboolean = .TRUE. 2253 if (fdf_output) write(fdf_out,'(l10,5x,a)') fdf_bboolean, & 2254 '# block designation by itself' 2255 endif 2256 2257 RETURN 2258 2259 CONTAINS 2260 2261 logical function is_true(valstr) result(a) 2262 character(len=*), intent(in) :: valstr 2263 a = leqi(valstr, 'yes') .or. leqi(valstr, 'true') .or. & 2264 leqi(valstr, '.true.') .or. leqi(valstr, 't') .or. & 2265 leqi(valstr, 'y') 2266 end function is_true 2267 2268 logical function is_false(valstr) result(a) 2269 character(len=*), intent(in) :: valstr 2270 a = leqi(valstr, 'no') .or. leqi(valstr, 'false') .or. & 2271 leqi(valstr, '.false.') .or. leqi(valstr, 'f') .or. & 2272 leqi(valstr, 'n') 2273 end function is_false 2274 2275!--------------------------------------------------------------------------- END 2276 END FUNCTION fdf_bboolean 2277 2278! 2279! Returns a single precision value associated with label 'label', 2280! or the default value if label is not found in the fdf file. 2281! Optionally can return a pointer to the line found. 2282! Note that integers on the line are also accepted 2283! 2284 FUNCTION fdf_single(label, default, line) 2285 implicit none 2286!--------------------------------------------------------------- Input Variables 2287 character(*) :: label 2288 real(sp) :: default 2289 2290!-------------------------------------------------------------- Output Variables 2291 real(sp) :: fdf_single 2292 type(line_dlist), pointer, optional :: line 2293 2294!--------------------------------------------------------------- Local Variables 2295 character(80) :: msg 2296 type(line_dlist), pointer :: mark 2297 2298!------------------------------------------------------------------------- BEGIN 2299! Prevents using FDF routines without initialize 2300 if (.not. fdf_started) then 2301 call die('FDF module: fdf_single', 'FDF subsystem not initialized', & 2302 THIS_FILE, __LINE__, fdf_err) 2303 endif 2304 2305 if (fdf_locate(label, mark)) then 2306 if (.not. match(mark%pline, 'lv')) then 2307 write(msg,*) 'no real value for ', label 2308 call die('FDF module: fdf_single', msg, THIS_FILE, __LINE__, fdf_err) 2309 endif 2310 fdf_single = values(mark%pline, 1, 1) 2311 if (fdf_output) write(fdf_out,'(a,5x,g20.10)') label, fdf_single 2312 else 2313 fdf_single = default 2314 if (fdf_output) write(fdf_out,'(a,5x,g20.10,5x,a)') label, default, '# default value' 2315 endif 2316 2317 if (PRESENT(line)) line = mark 2318 2319 RETURN 2320!--------------------------------------------------------------------------- END 2321 END FUNCTION fdf_single 2322 2323! 2324! Returns a double precision value associated with label 'label', 2325! or the default value if label is not found in the fdf file. 2326! Optionally can return a pointer to the line found. 2327! Note that integers on the line are also accepted 2328! 2329 FUNCTION fdf_double(label, default, line) 2330 implicit none 2331!--------------------------------------------------------------- Input Variables 2332 character(*) :: label 2333 real(dp) :: default 2334 2335!-------------------------------------------------------------- Output Variables 2336 real(dp) :: fdf_double 2337 type(line_dlist), pointer, optional :: line 2338 2339!--------------------------------------------------------------- Local Variables 2340 character(80) :: msg 2341 type(line_dlist), pointer :: mark 2342 2343!------------------------------------------------------------------------- BEGIN 2344! Prevents using FDF routines without initialize 2345 if (.not. fdf_started) then 2346 call die('FDF module: fdf_double', 'FDF subsystem not initialized', & 2347 THIS_FILE, __LINE__, fdf_err) 2348 endif 2349 2350 if (fdf_locate(label, mark)) then 2351 if (.not. match(mark%pline, 'lv')) then 2352 write(msg,*) 'no real value for ', label 2353 call die('FDF module: fdf_double', msg, THIS_FILE, __LINE__, fdf_err) 2354 endif 2355 fdf_double = values(mark%pline, 1, 1) 2356 if (fdf_output) write(fdf_out,'(a,5x,g20.10)') label, fdf_double 2357 else 2358 fdf_double = default 2359 if (fdf_output) write(fdf_out,'(a,5x,g20.10,5x,a)') label, default, '# default value' 2360 endif 2361 2362 if (PRESENT(line)) line = mark 2363 2364 RETURN 2365!--------------------------------------------------------------------------- END 2366 END FUNCTION fdf_double 2367 2368! 2369! Returns a double precision value associated with label 'label', 2370! or the default value if label is not found in the fdf file. 2371! Converts the units to defunit. 2372! Optionally can return a pointer to the line found. 2373! 2374 FUNCTION fdf_physical(label, default, defunit, line) 2375 implicit none 2376!--------------------------------------------------------------- Input Variables 2377 character(*) :: label, defunit 2378 real(dp) :: default 2379 2380!-------------------------------------------------------------- Output Variables 2381 real(dp) :: fdf_physical 2382 type(line_dlist), pointer, optional :: line 2383 2384!--------------------------------------------------------------- Local Variables 2385 character(10) :: unitstr 2386 character(80) :: msg 2387 real(dp) :: value 2388 type(line_dlist), pointer :: mark 2389 2390!------------------------------------------------------------------------- BEGIN 2391! Prevents using FDF routines without initialize 2392 if (.not. fdf_started) then 2393 call die('FDF module: fdf_physical', 'FDF subsystem not initialized', & 2394 THIS_FILE, __LINE__, fdf_err) 2395 endif 2396 2397! Label found 2398 if (fdf_locate(label, mark)) then 2399 if (.not. match(mark%pline, 'lv')) then 2400 write(msg,*) 'no real value for ', label 2401 call die('FDF module: fdf_physical', msg, THIS_FILE, __LINE__, fdf_err) 2402 endif 2403 2404! Label with value 2405 value = values(mark%pline, 1, 1) 2406 fdf_physical = value 2407 2408! Look for unit 2409 if (.not. match(mark%pline, 'lvn')) then 2410 write(msg,*) 'no unit specified for ', label 2411 call die('FDF module: fdf_physical', msg, THIS_FILE, __LINE__, fdf_err) 2412 endif 2413 2414 unitstr = names(mark%pline, 1, 2) 2415 if (.not. leqi(unitstr, defunit)) & 2416 fdf_physical = value * fdf_convfac(unitstr, defunit) 2417 2418 if (fdf_output) write(fdf_out,'(a,5x,g20.10,1x,a10)') label, fdf_physical, defunit 2419 if (fdf_output) write(fdf_out,'(a,a,5x,g20.10,1x,a10)') & 2420 '# above item originally: ', label, value, unitstr 2421 else 2422 fdf_physical = default 2423 if (fdf_output) write(fdf_out,'(a,5x,g20.10,1x,a,5x,a)') & 2424 label, default, defunit, '# default value' 2425 endif 2426 2427 if (PRESENT(line)) line = mark 2428 2429 RETURN 2430!--------------------------------------------------------------------------- END 2431 END FUNCTION fdf_physical 2432 2433! 2434! Returns a double precision value from a block-line after a certain input value 2435! or the default value if label is not found in the fdf file. 2436! Converts the units to defunit. 2437! 2438 FUNCTION fdf_bphysical(pline, default, defunit, after) 2439 implicit none 2440!--------------------------------------------------------------- Input Variables 2441 type(parsed_line), pointer :: pline 2442 real(dp) :: default 2443 character(*) :: defunit 2444 integer(ip), intent(in), optional :: after 2445 2446!-------------------------------------------------------------- Output Variables 2447 real(dp) :: fdf_bphysical 2448 2449!--------------------------------------------------------------- Local Variables 2450 character(10) :: unitstr 2451 character(80) :: msg 2452 real(dp) :: value 2453 type(line_dlist), pointer :: mark 2454 2455!------------------------------------------------------------------------- BEGIN 2456! Prevents using FDF routines without initialize 2457 if (.not. fdf_started) then 2458 call die('FDF module: fdf_bphysical', 'FDF subsystem not initialized', & 2459 THIS_FILE, __LINE__, fdf_err) 2460 endif 2461 2462 if (.not. match(pline, 'vn', after)) then 2463 write(msg,*) 'no real value for line: '//pline%line 2464 call die('FDF module: fdf_bphysical', msg, THIS_FILE, & 2465 __LINE__, fdf_err) 2466 endif 2467 2468 ! get value in block-line 2469 value = values(pline, 1, after) 2470 2471 ! get unit in block-line 2472 unitstr = names(pline, 1, after) 2473 if ( leqi(unitstr, defunit) ) then 2474 fdf_bphysical = value 2475 else 2476 fdf_bphysical = value * fdf_convfac(unitstr, defunit) 2477 end if 2478 2479 if ( fdf_output ) then 2480 if ( present(after) ) then 2481 write(fdf_out,'(5x,g20.10,1x,a10,1x,i0)') fdf_bphysical, & 2482 defunit, after 2483 else 2484 write(fdf_out,'(5x,g20.10,1x,a10)') fdf_bphysical, defunit 2485 end if 2486 write(fdf_out,'(a,a,5x,g20.10,1x,a10)') & 2487 '# above item on line: ', pline%line 2488 end if 2489 2490!--------------------------------------------------------------------------- END 2491 END FUNCTION fdf_bphysical 2492 2493! 2494! Returns conversion factor between a subset of physical units 2495! Written by j.m.soler. dec'96. 2496! Modified by alberto garcia, jan'97. 2497! Added more units by Nick Papior, Aug'17. 2498! 2499 FUNCTION fdf_convfac(from, to) 2500 implicit none 2501!--------------------------------------------------------------- Input Variables 2502 character(*) :: from, to 2503 2504!-------------------------------------------------------------- Output Variables 2505 real(dp) :: fdf_convfac 2506 2507!--------------------------------------------------------------- Local Variables 2508 character(80) :: msg 2509 integer(ip) :: iu, ifrom, ito 2510 2511!------------------------------------------------------------------------- BEGIN 2512 2513! 2514! We allow case variations in the units. this could be dangerous 2515! (meV --> MeV!!) in real life, but not in this restricted 2516! field. 2517! 2518 integer(ip), parameter :: nu = 78 2519 character(8) :: dimm(nu) 2520 character(10) :: name(nu) 2521 real(dp) :: unit(nu) 2522 data (dimm(iu), name(iu), unit(iu), iu=1, 3) / & 2523 'mass ', 'g ', 1.d-3, & 2524 'mass ', 'kg ', 1.d0, & 2525 'mass ', 'amu ', 1.66054d-27 / 2526 2527 data (dimm(iu), name(iu), unit(iu), iu=4, 9) / & 2528 'length ', 'm ', 1.d0, & 2529 'length ', 'cm ', 1.d-2, & 2530 'length ', 'nm ', 1.d-9, & 2531 'length ', 'pm ', 1.d-12, & 2532 'length ', 'ang ', 1.d-10, & 2533 'length ', 'bohr ', 0.529177d-10 / 2534 2535 data (dimm(iu), name(iu), unit(iu), iu=10, 19) / & 2536 'energy ', 'j ', 1.d0, & 2537 'energy ', 'kj ', 1.d3, & 2538 'energy ', 'erg ', 1.d-7, & 2539 'energy ', 'mev ', 1.60219d-22, & 2540 'energy ', 'ev ', 1.60219d-19, & 2541 'energy ', 'mry ', 2.17991d-21, & 2542 'energy ', 'ry ', 2.17991d-18, & 2543 'energy ', 'mha ', 4.35982d-21, & 2544 'energy ', 'mhartree ', 4.35982d-21, & 2545 'energy ', 'ha ', 4.35982d-18 / 2546 data (dimm(iu), name(iu), unit(iu), iu=20, 29) / & 2547 'energy ', 'hartree ', 4.35982d-18, & 2548 'energy ', 'k ', 1.38066d-23, & 2549 'energy ', 'kelvin ', 1.38066d-23, & 2550 'energy ', 'kcal/mol ', 6.94780d-21, & 2551 'energy ', 'kj/mol ', 1.6606d-21, & 2552 'energy ', 'hz ', 6.6262d-34, & 2553 'energy ', 'thz ', 6.6262d-22, & 2554 'energy ', 'cm-1 ', 1.986d-23, & 2555 'energy ', 'cm^-1 ', 1.986d-23, & 2556 'energy ', 'cm**-1 ', 1.986d-23 / 2557 2558 data (dimm(iu), name(iu), unit(iu), iu=30, 39) / & 2559 'time ', 's ', 1.d0, & 2560 'time ', 'ns ', 1.d-9, & 2561 'time ', 'ps ', 1.d-12, & 2562 'time ', 'fs ', 1.d-15, & 2563 'time ', 'min ', 60.d0, & 2564 'time ', 'mins ', 60.d0, & 2565 'time ', 'hour ', 3600.d0, & 2566 'time ', 'hours ', 3600.d0, & 2567 'time ', 'day ', 86400.d0, & 2568 'time ', 'days ', 86400.d0 / 2569 2570 data (dimm(iu), name(iu), unit(iu), iu=40, 43) / & 2571 'force ', 'n ', 1.d0, & 2572 'force ', 'ev/ang ', 1.60219d-9, & 2573 'force ', 'ry/bohr ', 4.11943d-8, & 2574 'force ', 'ha/bohr ', 8.23886d-08 / 2575 2576 data (dimm(iu), name(iu), unit(iu), iu=44, 52) / & 2577 'pressure', 'pa ', 1.d0, & 2578 'pressure', 'gpa ', 1.d9, & 2579 'pressure', 'atm ', 1.01325d5, & 2580 'pressure', 'bar ', 1.d5, & 2581 'pressure', 'mbar ', 1.d11, & 2582 'pressure', 'ev/ang**3 ', 1.60219d11, & 2583 'pressure', 'ev/ang^3 ', 1.60219d11, & 2584 'pressure', 'ry/bohr**3', 1.47108d13, & 2585 'pressure', 'ry/bohr^3 ', 1.47108d13 / 2586 2587 data (dimm(iu), name(iu), unit(iu), iu=53, 54) / & 2588 'charge ', 'c ', 1.d0, & 2589 'charge ', 'e ', 1.602177d-19 / 2590 2591 data (dimm(iu), name(iu), unit(iu), iu=55, 59) / & 2592 'dipole ', 'c*m ', 1.d0, & 2593 'dipole ', 'd ', 3.33564d-30, & 2594 'dipole ', 'debye ', 3.33564d-30, & 2595 'dipole ', 'e*bohr ', 8.47835d-30, & 2596 'dipole ', 'e*ang ', 1.602177d-29 / 2597 2598 data (dimm(iu), name(iu), unit(iu), iu=60, 61) / & 2599 'mominert', 'kg*m**2 ', 1.d0, & 2600 'mominert', 'ry*fs**2 ', 2.17991d-48 / 2601 2602 data (dimm(iu), name(iu), unit(iu), iu=62, 68) / & 2603 'efield ', 'v/m ', 1.d0, & 2604 'efield ', 'v/nm ', 1.d9, & 2605 'efield ', 'v/ang ', 1.d10, & 2606 'efield ', 'v/bohr ', 1.8897268d10, & 2607 'efield ', 'ry/bohr/e ', 2.5711273d11, & 2608 'efield ', 'ha/bohr/e ', 5.1422546d11, & 2609 'efield ', 'har/bohr/e', 5.1422546d11 / 2610 2611 data (dimm(iu), name(iu), unit(iu), iu=69, 70) / & 2612 'angle ', 'deg ', 1.d0, & 2613 'angle ', 'rad ', 5.72957795d1 / 2614 2615 data (dimm(iu), name(iu), unit(iu), iu=71, 78) / & 2616 'torque ', 'mev/deg ', 1.0d-3, & 2617 'torque ', 'mev/rad ', 1.745533d-5, & 2618 'torque ', 'ev/deg ', 1.0d0, & 2619 'torque ', 'ev/rad ', 1.745533d-2, & 2620 'torque ', 'mry/deg ', 13.6058d-3, & 2621 'torque ', 'mry/rad ', 0.237466d-3, & 2622 'torque ', 'ry/deg ', 13.6058d0, & 2623 'torque ', 'ry/rad ', 0.237466d0 / 2624 2625! 2626 ifrom = 0 2627 ito = 0 2628 do iu= 1, nu 2629 if (leqi(name(iu), from)) ifrom = iu 2630 if (leqi(name(iu), to)) ito = iu 2631 end do 2632 2633 if (ifrom .eq. 0) then 2634 write(msg,*) 'unknown unit = ', from 2635 call die('FDF module: fdf_convfac', msg, THIS_FILE, __LINE__, fdf_err) 2636 endif 2637 2638 if (ito .eq. 0) then 2639 write(msg,*) 'unknown unit = ', to 2640 call die('FDF module: fdf_convfac', msg, THIS_FILE, __LINE__, fdf_err) 2641 endif 2642 2643 if (leqi(dimm(ifrom), dimm(ito))) then 2644 fdf_convfac = unit(ifrom) / unit(ito) 2645 else 2646 write(msg,*) 'unit''s physical dimensions don''t match: ', & 2647 from, ', ', to 2648 call die('FDF module: fdf_convfac', msg, THIS_FILE, __LINE__, fdf_err) 2649 endif 2650!--------------------------------------------------------------------------- END 2651 END FUNCTION fdf_convfac 2652 2653! 2654! Searches for label in the fdf hierarchy. If it appears the function 2655! returns .TRUE. and leaves mark pointer positioned at the line. 2656! Otherwise, it returns .FALSE. and mark points to NULL. 2657! 2658 FUNCTION fdf_locate(label, mark) 2659 implicit none 2660!--------------------------------------------------------------- Input Variables 2661 character(*) :: label 2662 2663!-------------------------------------------------------------- Output Variables 2664 logical :: fdf_locate 2665 type(line_dlist), pointer :: mark 2666 2667!--------------------------------------------------------------- Local Variables 2668 character(80) :: strlabel 2669 2670!------------------------------------------------------------------------- BEGIN 2671 fdf_locate = .FALSE. 2672 2673! if (fdf_donothing) return 2674 2675 mark => file_in%first 2676 do while ((.not. fdf_locate) .and. (ASSOCIATED(mark))) 2677 2678 if (match(mark%pline, 'l')) then 2679 strlabel = labels(mark%pline) 2680 2681 if (labeleq(strlabel, label, fdf_log)) then 2682 fdf_locate = .TRUE. 2683 else 2684 mark => mark%next 2685 endif 2686 else 2687 mark => mark%next 2688 endif 2689 enddo 2690 2691 if (.not. fdf_locate) NULLIFY(mark) 2692 2693 RETURN 2694!--------------------------------------------------------------------------- END 2695 END FUNCTION fdf_locate 2696 2697! 2698! Returns true or false whether or not the label 'label' is 2699! a block. 2700! I.e. it returns true if the line has the form bl, if not found, or not bl 2701! it returns false. 2702! 2703 FUNCTION fdf_isblock(label) 2704 implicit none 2705!--------------------------------------------------------------- Input Variables 2706 character(*) :: label 2707 2708!-------------------------------------------------------------- Output Variables 2709 logical :: fdf_isblock 2710 2711!--------------------------------------------------------------- Local Variables 2712 type(line_dlist), pointer :: mark 2713 character(80) :: strlabel 2714 2715!------------------------------------------------------------------------- BEGIN 2716! Prevents using FDF routines without initialize 2717 if (.not. fdf_started) then 2718 call die('FDF module: fdf_isblock', 'FDF subsystem not initialized', & 2719 THIS_FILE, __LINE__, fdf_err) 2720 endif 2721 2722 fdf_isblock = .false. 2723 2724 mark => file_in%first 2725 do while ( associated(mark) ) 2726 2727!!$ if ( match(mark%pline, 'l') ) then 2728!!$ strlabel = labels(mark%pline) 2729!!$ 2730!!$ if ( labeleq(strlabel, label, fdf_log) ) then 2731!!$ ! fdf has first-encounter acceptance. 2732!!$ ! I.e. for an input 2733!!$ ! Label_Name 1 2734!!$ ! %block Label_Name 2735!!$ ! 1 2736!!$ ! %endblock Label_Name 2737!!$ ! the former will be accepted first. 2738!!$ exit 2739!!$ end if 2740 2741 if ( match(mark%pline, 'bl') ) then 2742 strlabel = blocks(mark%pline) 2743 2744 if ( labeleq(strlabel, label, fdf_log) ) then 2745 fdf_isblock = .true. 2746 exit 2747 end if 2748 end if 2749 2750 mark => mark%next 2751 end do 2752 2753 if (fdf_output) write(fdf_out,'(a,5x,l10)') "#:block? " // label, fdf_isblock 2754 2755 RETURN 2756!--------------------------------------------------------------------------- END 2757 END FUNCTION fdf_isblock 2758 2759! 2760! Searches for block label in the fdf hierarchy. If it appears returns 2761! .TRUE. and leaves block mark pointer positioned at the first line. 2762! Otherwise, it returns .FALSE. and block mark points to NULL. 2763! 2764 FUNCTION fdf_block(label, bfdf) 2765 implicit none 2766!--------------------------------------------------------------- Input Variables 2767 character(*) :: label 2768 2769!-------------------------------------------------------------- Output Variables 2770 logical :: fdf_block 2771 type(block_fdf) :: bfdf 2772 2773!--------------------------------------------------------------- Local Variables 2774 character(80) :: strlabel 2775 2776!------------------------------------------------------------------------- BEGIN 2777! Prevents using FDF routines without initialize 2778 if (.not. fdf_started) then 2779 call die('FDF module: fdf_block', 'FDF subsystem not initialized', & 2780 THIS_FILE, __LINE__, fdf_err) 2781 endif 2782 2783 fdf_block = .FALSE. 2784 2785 bfdf%mark => file_in%first 2786 do while ((.not. fdf_block) .and. (ASSOCIATED(bfdf%mark))) 2787 2788 if (match(bfdf%mark%pline, 'bl')) then 2789 strlabel = blocks(bfdf%mark%pline) 2790 2791 if (labeleq(strlabel, label, fdf_log)) then 2792 fdf_block = .TRUE. 2793 bfdf%label = label 2794 2795 if (fdf_output) write(fdf_out,'(a,a)') '%block ', TRIM(label) 2796 endif 2797 endif 2798 2799 bfdf%mark => bfdf%mark%next 2800 enddo 2801 2802 if (.not. fdf_block) NULLIFY(bfdf%mark) 2803 2804 RETURN 2805!--------------------------------------------------------------------------- END 2806 END FUNCTION fdf_block 2807 2808! 2809! Get successive parsed lines from block returning 2810! .TRUE. while more lines exist in the block bfdf. 2811! 2812 FUNCTION fdf_bline(bfdf, pline) 2813 implicit none 2814!--------------------------------------------------------------- Input Variables 2815 type(block_fdf) :: bfdf 2816 2817!-------------------------------------------------------------- Output Variables 2818 logical :: fdf_bline 2819 type(parsed_line), pointer :: pline 2820 2821!--------------------------------------------------------------- Local Variables 2822 character(80) :: strlabel 2823 2824!------------------------------------------------------------------------- BEGIN 2825! Prevents using FDF routines without initialize 2826 if (.not. fdf_started) then 2827 call die('FDF module: fdf_bline', 'FDF subsystem not initialized', & 2828 THIS_FILE, __LINE__, fdf_err) 2829 endif 2830 2831 if (.not. ASSOCIATED(bfdf%mark)) then 2832 call die('FDF module: fdf_bline', 'block_fdf structure not initialized', & 2833 THIS_FILE, __LINE__, fdf_err) 2834 endif 2835 2836 fdf_bline = .TRUE. 2837 2838! If we are in the head of the block move to the content 2839 if (match(bfdf%mark%pline, 'bl')) then 2840 strlabel = blocks(bfdf%mark%pline) 2841 2842 if (labeleq(strlabel, bfdf%label, fdf_log)) then 2843 bfdf%mark => bfdf%mark%next 2844 2845 if (fdf_output) write(fdf_out,'(a,a)') '%block ', TRIM(bfdf%label) 2846 endif 2847 endif 2848 2849 if (match(bfdf%mark%pline, 'el')) then 2850 strlabel = endblocks(bfdf%mark%pline) 2851 2852 if (labeleq(strlabel, bfdf%label, fdf_log)) then 2853 fdf_bline = .FALSE. 2854 NULLIFY(pline) 2855 2856 if (fdf_output) write(fdf_out,'(a,a)') '%endblock ', TRIM(bfdf%label) 2857 endif 2858 endif 2859 2860 if (fdf_bline) then 2861 if (fdf_output) write(fdf_out,'(1x,a)') TRIM(bfdf%mark%str) 2862 2863 pline => bfdf%mark%pline 2864 bfdf%mark => bfdf%mark%next 2865 endif 2866 2867 RETURN 2868!--------------------------------------------------------------------------- END 2869 END FUNCTION fdf_bline 2870 2871 2872 2873! 2874! Backspace to the previous physical line in the block 2875! returning .TRUE. while more lines exist in the block bfdf. 2876! 2877 FUNCTION fdf_bbackspace(bfdf,pline) 2878 implicit none 2879!--------------------------------------------------------------- Input Variables 2880 type(block_fdf) :: bfdf 2881 2882!-------------------------------------------------------------- Output Variables 2883 logical :: fdf_bbackspace 2884 type(parsed_line), pointer, optional :: pline 2885 2886!--------------------------------------------------------------- Local Variables 2887 character(80) :: strlabel 2888 2889!------------------------------------------------------------------------- BEGIN 2890! Prevents using FDF routines without initialize 2891 if (.not. fdf_started) then 2892 call die('FDF module: fdf_bbackspace', 'FDF subsystem not initialized', & 2893 THIS_FILE, __LINE__, fdf_err) 2894 endif 2895 2896 if (.not. ASSOCIATED(bfdf%mark)) then 2897 call die('FDF module: fdf_bbackspace', 'block_fdf structure not initialized', & 2898 THIS_FILE, __LINE__, fdf_err) 2899 endif 2900 2901 fdf_bbackspace = .TRUE. 2902 2903! If we are in the bottom of the block move to the content 2904 2905 if (match(bfdf%mark%pline, 'el')) then 2906 2907 strlabel = endblocks(bfdf%mark%pline) 2908 2909 if (labeleq(strlabel, bfdf%label, fdf_log)) then 2910 bfdf%mark => bfdf%mark%prev 2911 2912 if (fdf_output) write(fdf_out,'(1x,a)') "#:(Backspace to) " // "|" // & 2913 TRIM(bfdf%mark%str) // "|" 2914 endif 2915 2916! If we are at the head we cannot backspace 2917 2918 else if (match(bfdf%mark%pline, 'bl')) then 2919 strlabel = blocks(bfdf%mark%pline) 2920 2921 if (labeleq(strlabel, bfdf%label, fdf_log)) then 2922 fdf_bbackspace = .FALSE. 2923 if (fdf_output) write(fdf_out,'(1x,a)') "#:(Cannot backspace) " // "|" // & 2924 TRIM(bfdf%mark%str) // "|" 2925 endif 2926 2927 else 2928 2929 bfdf%mark => bfdf%mark%prev 2930 if (fdf_output) write(fdf_out,'(1x,a)') "#:(Backspace to) " // "|" // & 2931 TRIM(bfdf%mark%str) // "|" 2932 endif 2933 2934 if ( present(pline) ) pline => bfdf%mark%pline 2935 2936 RETURN 2937!--------------------------------------------------------------------------- END 2938 END FUNCTION fdf_bbackspace 2939 2940! 2941! Moves the pointer of the working line (bfdf%mark) 2942! to the beginning of the block 'label' structure. 2943! 2944 SUBROUTINE fdf_brewind(bfdf) 2945 implicit none 2946!-------------------------------------------------------------- Output Variables 2947 type(block_fdf) :: bfdf 2948 2949!--------------------------------------------------------------- Local Variables 2950 character(80) :: msg 2951 2952!------------------------------------------------------------------------- BEGIN 2953! Prevents using FDF routines without initialize 2954 if (.not. fdf_started) then 2955 call die('FDF module: fdf_brewind', 'FDF subsystem not initialized', & 2956 THIS_FILE, __LINE__, fdf_err) 2957 endif 2958 2959 if (.not. ASSOCIATED(bfdf%mark)) then 2960 call die('FDF module: fdf_brewind', 'block_fdf structure not initialized', & 2961 THIS_FILE, __LINE__, fdf_err) 2962 endif 2963 2964 if (.not. fdf_block(bfdf%label, bfdf)) then 2965 write(msg,*) 'Block ', bfdf%label, ' not found in FDF structure' 2966 call die('FDF module: fdf_brewind', msg, & 2967 THIS_FILE, __LINE__, fdf_err) 2968 endif 2969 2970 RETURN 2971!--------------------------------------------------------------------------- END 2972 END SUBROUTINE fdf_brewind 2973 2974! 2975! Closes the opened block by looping the remaining lines of the working line. 2976! This is only needed to complete the fdf-*.log files output for direct 2977! usage later. It does nothing internally. 2978! 2979 SUBROUTINE fdf_bclose(bfdf) 2980 implicit none 2981!-------------------------------------------------------------- Output Variables 2982 type(block_fdf) :: bfdf 2983 2984!--------------------------------------------------------------- Local Variables 2985 type(parsed_line), pointer :: pline 2986 integer(ip) :: i 2987 character(80) :: msg 2988 2989!------------------------------------------------------------------------- BEGIN 2990! Prevents using FDF routines without initialize 2991 if (.not. fdf_started) then 2992 call die('FDF module: fdf_bclose', 'FDF subsystem not initialized', & 2993 THIS_FILE, __LINE__, fdf_err) 2994 endif 2995 2996 ! Quick return (no need for errors) 2997 if ( .not. associated(bfdf%mark) ) return 2998 2999 ! This should hopefully discourage compilers to optimize the loop away... 3000 i = 0 3001 do while ( fdf_bline(bfdf, pline) ) 3002 i = i + fdf_bnvalues(pline) 3003 end do 3004 write(msg,'(a,i10)') 'Block ', i 3005 3006 RETURN 3007!--------------------------------------------------------------------------- END 3008 END SUBROUTINE fdf_bclose 3009 3010 3011! 3012! Count number of lines with an optional specification. 3013! I.e. this will read through the block and return the number of lines in the 3014! block which matches the morphology (morph) 3015! This may be used to easily digest number of non-empty lines in the block. 3016! Note that a match on the morphology only compares the number of ID's in 3017! morph. I.e. a line with 'vvvil' will match 'vvvi'. 3018! 3019 FUNCTION fdf_block_linecount(label, morph) 3020 implicit none 3021!--------------------------------------------------------------- Input Variables 3022 character(len=*) :: label 3023 character(len=*), optional :: morph 3024!-------------------------------------------------------------- Output Variables 3025 integer(ip) :: fdf_block_linecount 3026 3027!--------------------------------------------------------------- Local Variables 3028 type(block_fdf) :: bfdf 3029 type(parsed_line), pointer :: pline 3030 logical :: orig_fdf_output 3031 3032!------------------------------------------------------------------------- BEGIN 3033! Prevents using FDF routines without initialize 3034 if (.not. fdf_started) then 3035 call die('FDF module: fdf_block_linecount', 'FDF subsystem not initialized', & 3036 THIS_FILE, __LINE__, fdf_err) 3037 endif 3038 3039 ! Store the fdf_output variable (suppress writing to log) 3040 orig_fdf_output = fdf_output 3041 fdf_output = .false. 3042 3043 ! Find the block and search for morhp 3044 fdf_block_linecount = 0 3045 if ( fdf_block(label, bfdf) ) then 3046 3047 do while ( fdf_bline(bfdf, pline) ) 3048 if ( present(morph) ) then 3049 if ( fdf_bmatch(pline, morph) ) then 3050 fdf_block_linecount = fdf_block_linecount + 1 3051 end if 3052 else 3053 fdf_block_linecount = fdf_block_linecount + 1 3054 end if 3055 end do 3056 3057 end if 3058 3059 ! Restore output 3060 fdf_output = orig_fdf_output 3061 3062 if ( fdf_output ) then 3063 if ( present(morph) ) then 3064 write(fdf_out,'(3a,3x,i0)') "#:block-line-count? ", & 3065 trim(label), ' ('//trim(morph)//')', fdf_block_linecount 3066 else 3067 write(fdf_out,'(2a,3x,i0)') "#:block-line-count? ", & 3068 trim(label), fdf_block_linecount 3069 end if 3070 end if 3071 3072 RETURN 3073!--------------------------------------------------------------------------- END 3074 END FUNCTION fdf_block_linecount 3075 3076! 3077! Check if label is defined 3078! 3079 logical FUNCTION fdf_defined(label) 3080 implicit none 3081!--------------------------------------------------------------- Input Variables 3082 character(*) :: label 3083 3084!--------------------------------------------------------------- Local Variables 3085 type(line_dlist), pointer :: mark 3086 3087!--------------------------------------------------------------------- BEGIN 3088 ! First, check whether a single label exists: 3089 fdf_defined = fdf_locate(label, mark) 3090 if (.not. fdf_defined) then 3091 ! Check whether there is a block with that label 3092 fdf_defined = fdf_isblock(label) 3093 endif 3094 if ( fdf_output ) then 3095 write(fdf_out,'(a,5x,l10)') '#:defined? ' // label, fdf_defined 3096 endif 3097 3098 RETURN 3099!----------------------------------------------------------------------- END 3100 END FUNCTION fdf_defined 3101 3102! 3103! Output levels: 3104! level <= 0: nothing 3105! level = 1: standard 3106! 3107 SUBROUTINE fdf_setoutput(level,fileout_in) 3108 implicit none 3109!------------------------------------------------------------- Input Variables 3110 integer(ip) :: level 3111 character(len=*), intent(in) :: fileout_in 3112 3113 3114 character(len=256) :: fileout 3115 3116 fileout = fileout_in 3117 if (level .le. 0) then 3118 if (fdf_output) then 3119 call io_close(fdf_out) 3120 fdf_output = .FALSE. 3121 endif 3122 else 3123 if (.not. fdf_output) then 3124#ifdef _MPI_ 3125 if (rank /= 0) then 3126 if (level .ge. 2) then 3127 fileout = "/tmp/" // trim(fileout) // "." // i2s(rank) 3128 call io_assign(fdf_out) 3129 open(fdf_out, file=fileout, form='formatted', & 3130 status='unknown') 3131 REWIND(fdf_out) 3132 fdf_output = .TRUE. 3133 endif 3134 else 3135 call io_assign(fdf_out) 3136 open(fdf_out, file=fileout, form='formatted', & 3137 status='unknown') 3138 REWIND(fdf_out) 3139 fdf_output = .TRUE. 3140 endif 3141#else 3142 call io_assign(fdf_out) 3143 open(fdf_out, file=fileout, form='formatted', & 3144 status='unknown') 3145 REWIND(fdf_out) 3146 fdf_output = .TRUE. 3147#endif 3148 endif 3149 endif 3150!----------------------------------------------------------------------- END 3151 END SUBROUTINE fdf_setoutput 3152 3153! 3154! Debugging levels: 3155! level <= 0: nothing 3156! level = 1: standard 3157! level >= 2: exhaustive 3158 3159 SUBROUTINE fdf_setdebug(level,filedebug) 3160 implicit none 3161!------------------------------------------------------------- Input Variables 3162 integer(ip) :: level 3163 character(len=*) :: filedebug 3164 3165!----------------------------------------------------------------------- BEGIN 3166 if (level .le. 0) then 3167 if (fdf_debug) then 3168 call io_close(fdf_log) 3169 fdf_debug = .FALSE. 3170 endif 3171 else 3172 if (.not. fdf_debug) then 3173 call io_assign(fdf_log) 3174#ifdef _MPI_ 3175 if (rank /= 0) then 3176 filedebug = "/tmp/" // trim(filedebug) // "." // i2s(rank) 3177 endif 3178#endif 3179 open(fdf_log, file=filedebug, form='formatted', & 3180 status='unknown') 3181 REWIND(fdf_log) 3182 fdf_debug = .TRUE. 3183 3184! Set logging/debugging info for PARSE module also 3185 call setlog(fdf_log) 3186 call setdebug(1) 3187 endif 3188 endif 3189 3190 fdf_debug2 = (level .ge. 2) 3191 3192 RETURN 3193!----------------------------------------------------------------------- END 3194 END SUBROUTINE fdf_setdebug 3195 3196! 3197! For handling deprecated labels. 3198! Also there is an optional "newlabel" if it has been changed into 3199! a new label. 3200! 3201 subroutine fdf_deprecated(label,newlabel) 3202 implicit none 3203!--------------------------------------------------------------- Input Variables 3204 character(*) :: label 3205 character(*) :: newlabel 3206 3207!------------------------------------------------------------------------- BEGIN 3208 if ( fdf_defined(label) ) then 3209 if (fdf_output) write(fdf_out,'(a)') "#**Warning: FDF symbol '"//trim(label)// & 3210 "' is deprecated." 3211 if ( fdf_defined(newlabel) ) then 3212 if (fdf_output) write(fdf_out,'(a)') "# FDF symbol '"//trim(newlabel)// & 3213 "' will be used instead." 3214 else 3215 if (fdf_output) write(fdf_out,'(a)') "# FDF symbol '"//trim(newlabel)// & 3216 "' replaces '"//trim(label)//"'." 3217 end if 3218 end if 3219 3220!--------------------------------------------------------------------------- END 3221 end subroutine fdf_deprecated 3222 3223! 3224! For handling obsoleted labels. 3225! 3226 subroutine fdf_obsolete(label) 3227 implicit none 3228!--------------------------------------------------------------- Input Variables 3229 character(*) :: label 3230 3231!------------------------------------------------------------------------- BEGIN 3232 if ( fdf_defined(label) ) then 3233 if (fdf_output) write(fdf_out,'(a)') "#**Warning: FDF symbol '"//trim(label)// & 3234 "' is obsolete." 3235 end if 3236 3237!--------------------------------------------------------------------------- END 3238 end subroutine fdf_obsolete 3239 3240 3241 subroutine serialize_fdf_struct(buffer) 3242 character, pointer :: buffer(:) 3243 3244 character(len=SERIALIZED_LENGTH) bufline 3245 type(line_dlist), pointer :: mark 3246 integer(ip) :: i, length, init, final 3247 3248 mark => file_in%first 3249 do i= 1, file_in%nlines 3250 call serialize_pline(mark%pline,bufline,length) 3251 init = (i-1)*SERIALIZED_LENGTH+1 3252 final = (i)*SERIALIZED_LENGTH 3253 call convert_string_to_array_of_chars(bufline,buffer(init:final)) 3254 mark => mark%next 3255 enddo 3256 end subroutine serialize_fdf_struct 3257 3258 subroutine recreate_fdf_struct(nlines,bufferFDF) 3259 character, pointer :: bufferFDF(:) 3260 integer, intent(in) :: nlines 3261 3262 character(len=SERIALIZED_LENGTH) bufline 3263 type(parsed_line), pointer :: pline 3264 integer(ip) :: i, init, final 3265 3266 do i= 1, nlines 3267 init = (i-1)*SERIALIZED_LENGTH+1 3268 final = (i)*SERIALIZED_LENGTH 3269 call convert_array_of_chars_to_string(bufferFDF(init:final),bufline) 3270 allocate(pline) 3271 call recreate_pline(pline,bufline) 3272 call fdf_addtoken(pline%line,pline) 3273 enddo 3274 ! print *, "Processed: ", file_in%nlines, " lines." 3275 3276 end subroutine recreate_fdf_struct 3277 3278!================================================================= 3279!--------- Obsolete routines ---------------- 3280! 3281! The owner node of the input file sends data file to the 3282! other processes in the MPI communicator. 3283! 3284#ifdef CLUSTER 3285 SUBROUTINE fdf_sendInput() 3286 use mpi_siesta 3287 implicit none 3288!--------------------------------------------------------------- Local Variables 3289 character, pointer :: bufferFDF(:) => null() 3290 integer(ip) :: i, j, k, ierr 3291 type(line_dlist), pointer :: mark 3292 3293!------------------------------------------------------------------------- BEGIN 3294 ALLOCATE(bufferFDF(file_in%nlines*MAX_LENGTH), stat=ierr) 3295 if (ierr .ne. 0) then 3296 call die('FDF module: fdf_sendInput', 'Error allocating bufferFDF', & 3297 THIS_FILE, __LINE__, fdf_err, rc=ierr) 3298 endif 3299 3300 mark => file_in%first 3301 do i= 1, file_in%nlines*MAX_LENGTH, MAX_LENGTH 3302 bufferFDF(i:i+MAX_LENGTH-1) = s2arr(mark%str) 3303 mark => mark%next 3304 enddo 3305 3306 call MPI_Bcast(file_in%nlines, 1, & 3307 MPI_INTEGER, rank, MPI_COMM_WORLD, ierr) 3308 3309 if (ierr .ne. MPI_SUCCESS) then 3310 call die('FDF module: fdf_sendInput', 'Error Broadcasting nlines.' // & 3311 'Terminating.', THIS_FILE, __LINE__, fdf_err, rc=ierr) 3312 endif 3313 3314 call MPI_Bcast(bufferFDF, file_in%nlines*MAX_LENGTH, & 3315 MPI_CHARACTER, rank, MPI_COMM_WORLD, ierr) 3316 if (ierr .ne. MPI_SUCCESS) then 3317 call die('FDF module: fdf_sendInput', 'Error Broadcasting bufferFDF.' // & 3318 'Terminating.', THIS_FILE, __LINE__, fdf_err, rc=ierr) 3319 endif 3320 3321 DEALLOCATE(bufferFDF) 3322 3323 RETURN 3324!--------------------------------------------------------------------------- END 3325 END SUBROUTINE fdf_sendInput 3326#endif 3327 3328! 3329! Remaining nodes recieve the input data file from the owner 3330! of the FDF file inside the cluster. 3331! 3332#ifdef CLUSTER 3333 SUBROUTINE fdf_recvInput(root, filein, fileinTmp) 3334 use mpi_siesta 3335 implicit none 3336!--------------------------------------------------------------- Input Variables 3337 character(*) :: filein, fileinTmp 3338 integer(ip) :: root 3339 3340!--------------------------------------------------------------- Local Variables 3341 integer(ip) :: i, j, lun, ierr, nlines 3342 character, pointer :: bufferFDF(:) => null() 3343 character(len=10) :: fmt 3344!------------------------------------------------------------------------- BEGIN 3345 call MPI_Bcast(nlines, 1, & 3346 MPI_INTEGER, root, MPI_COMM_WORLD, ierr) 3347 if (ierr .ne. MPI_SUCCESS) then 3348 call die('FDF module: fdf_recvInput', 'Error Broadcasting nlines.' // & 3349 'Terminating.', THIS_FILE, __LINE__, fdf_err, rc=ierr) 3350 endif 3351 3352 ALLOCATE(bufferFDF(nlines*MAX_LENGTH), stat=ierr) 3353 if (ierr .ne. 0) then 3354 call die('FDF module: fdf_recvInput', 'Error allocating bufferFDF', & 3355 THIS_FILE, __LINE__, fdf_err, rc=ierr) 3356 endif 3357 3358 call MPI_Bcast(bufferFDF, nlines*MAX_LENGTH, & 3359 MPI_CHARACTER, root, MPI_COMM_WORLD, ierr) 3360 if (ierr .ne. MPI_SUCCESS) then 3361 call die('FDF module: fdf_recvInput', 'Error Broadcasting bufferFDF.' // & 3362 'Terminating.', THIS_FILE, __LINE__, fdf_err, rc=ierr) 3363 endif 3364 3365 call io_assign(lun) 3366 fileinTmp = TRIM(filein) // '.' // i2s(rank) 3367 open(unit=lun, file=fileinTmp, form='formatted', & 3368 status='unknown') 3369 3370 if (MAX_LENGTH.lt.10) then 3371 write(fmt,"(a1,I1,a2)") "(", MAX_LENGTH, "a)" 3372 else if (MAX_LENGTH.lt.100) then 3373 write(fmt,"(a1,I2,a2)") "(", MAX_LENGTH, "a)" 3374 else if (MAX_LENGTH.lt.1000) then 3375 write(fmt,"(a1,I3,a2)") "(", MAX_LENGTH, "a)" 3376 else 3377 write(fmt,"(a1,I4,a2)") "(", MAX_LENGTH, "a)" 3378 endif 3379 3380 do i= 1, nlines*MAX_LENGTH, MAX_LENGTH 3381 write(lun,fmt) (bufferFDF(j),j=i,i+MAX_LENGTH-1) 3382 enddo 3383 call io_close(lun) 3384 3385 DEALLOCATE(bufferFDF) 3386 3387 RETURN 3388!--------------------------------------------------------------------------- END 3389 END SUBROUTINE fdf_recvInput 3390#endif 3391END MODULE fdf 3392