TABLE OF CONTENTS


ABINIT/intagm [ Functions ]

[ Top ] [ Functions ]

NAME

 intagm

FUNCTION

 Search input 'string' for specific 'token'. Search depends on
 input dataset through 'jdtset'. Then, return the information
 mentioned after 'token'.
 See the "notes" section

COPYRIGHT

 Copyright (C) 1998-2017 ABINIT group (DCA, XG, GMR).
 This file is distributed under the terms of the
 GNU General Public License, see ~abinit/COPYING
 or http://www.gnu.org/copyleft/gpl.txt .
 For the initials of contributors, see ~abinit/doc/developers/contributors.txt .

INPUTS

  jdtset=see the notes section
  marr=dimension of the intarr and dprarr arrays, as declared in the calling subroutine.
  narr=actual size of array to be read in.
  string=character string containing 'tags' and data.
  token=character string for 'tag'.
  typevarphys= variable type (might indicate the physical meaning of for dimensionality purposes)
   'INT'=>integer
   'DPR'=>real(dp) (no special treatment)
   'LEN'=>real(dp) (expect a "length", identify bohr, au or angstrom,
       and return in au -atomic units=bohr- )
   'ENE'=>real(dp) (expect a "energy", identify Ha, hartree, eV, Ry, Rydberg)
   'LOG'=>integer, but read logical variable T,F,.true., or .false.
   'KEY'=>character, returned in key_value

OUTPUT

  intarr(1:narr), dprarr(1:narr)
   integer or real(dp) arrays, respectively (see typevarphys),
   into which data is read if typevarphys/='KEY'. Use these arrays even for scalars.
  tread is an integer : tread = 0 => no data was read
                        tread = 1 => data was read
  ds_input is an optional integer flag:
           ds_input = 0 => value was found which is not specific to jdtset
           ds_input > 0 => value was found which is specific to jdtset
   one could add more information, eg whether a ? or a : was used, etc...
   [key_value]=Stores the value of key if typevarphys=="KEY"
               String of len fnlen

NOTES

 If jdtset==0:

  Search compressed 'string' for blank//'token'//blank and
  read input data beside 'token', to be read into appropriate variable.
  For this routine to find a given token, the token has to be preceded
  and followed by blanks--i.e. the first token should not start out as
  the first character in the input file.  This is checked in the calling
  subroutine 'input'. Calls inread which performs internal read from
  specified string.  Also calls upper which maps characters to all upper case.
  Also checks whether there is an occurence of blank//'token'//digit,
  in which case the input file might be erroneous, so stops.

 If jdtset is a positive number:

  (1) First search for modified string, blank//'token'//jdtset//blank

  (2a) if the occurence of (1) is not found,
       look for other modified strings,
       blank//'token'//'?'//unities//blank
       or
       blank//'token'//dozens//'?'//blank
       (issue an error message if more than one occurs)
       where jdtset=dozens*10+unities (decimal decomposition of jdtset)
       if one of them exists, just take the value
       Note that unities is a one-digit number, while dozens might be bigger than 9.

  (2b-2c) search for a series, with the following tokens :
       (issue an error message if more than one occurs, or
       goto (3) if none exist)

      blank//'token'//':'//blank
      if it exists, then a series might have been defined in the input file
      must thus find either the increment, blank//'token'//'+'//blank,
      or the multiplicative factor, blank//'token'//'*'//blank

      blank//'token'//'?'//':'//blank
      if it exists, then a series for the inner loop
      might have been defined in the input file
      must thus find either the increment, blank//'token'//'?'//'+'//blank,
      or the multiplicative factor, blank//'token'//'?'//'*'//blank

      blank//'token'//':'//'?'//blank
      if it exists, then a series for the outer loop
      might have been defined in the input file
      must thus find either the increment, blank//'token'//'+'//'?'//blank,
      or the multiplicative factor, blank//'token'//'*'//'?'//blank

  (3) if neither (1) nor (2) are found, search for the 'normal'
       string, blank//'token'//blank

PARENTS

      ingeo,ingeobld,inkpts,inqpt,invacuum,invars0,invars1,invars2
      m_ab7_invars_f90,m_anaddb_dataset,m_band2eps_dataset,m_ingeo_img
      m_multibinit_dataset,macroin,mpi_setup,parsefile,ujdet

CHILDREN

      appdig,inarray,inupper,wrtout

SOURCE

110 #if defined HAVE_CONFIG_H
111 #include "config.h"
112 #endif
113 
114 #include "abi_common.h"
115 
116 subroutine intagm(dprarr,intarr,jdtset,marr,narr,string,token,tread,typevarphys,ds_input,key_value)
117 
118  use defs_basis
119  use m_profiling_abi
120  use m_errors
121 
122  use m_fstrings,        only : sjoin, itoa
123 
124 !This section has been created automatically by the script Abilint (TD).
125 !Do not modify the following lines by hand.
126 #undef ABI_FUNC
127 #define ABI_FUNC 'intagm'
128  use interfaces_14_hidewrite
129  use interfaces_32_util
130  use interfaces_42_parser, except_this_one => intagm
131 !End of the abilint section
132 
133  implicit none
134 
135 !Arguments ------------------------------------
136 !scalars
137  integer,intent(in) :: jdtset,marr,narr
138  integer,intent(out) :: tread
139  integer,intent(out),optional :: ds_input
140  character(len=*),intent(in) :: string
141  character(len=*),intent(in) :: token
142  character(len=*),intent(in) :: typevarphys
143  character(len=fnlen),optional,intent(out) :: key_value
144 !arrays
145  integer,intent(inout) :: intarr(marr) !vz_i
146  real(dp),intent(inout) :: dprarr(marr) !vz_i
147 
148 !Local variables-------------------------------
149  character(len=1), parameter :: blank=' '
150 !scalars
151  integer :: b1,cs1len,cslen,dozens,ier,itoken,itoken1,itoken2,itoken2_1colon
152  integer :: itoken2_1plus,itoken2_1times,itoken2_2colon,itoken2_2plus
153  integer :: itoken2_2times,itoken2_colon,itoken2_plus,itoken2_times
154  integer :: itoken_1colon,itoken_1plus,itoken_1times,itoken_2colon,itoken_2plus
155  integer :: itoken_2times,itoken_colon,itoken_plus,itoken_times,number,opttoken
156  integer :: sum_token,toklen,trial_cslen,trial_jdtset,unities
157  integer :: ds_input_
158  character(len=4) :: appen
159  character(len=3) :: typevar
160  character(len=500) :: message
161  character(len=fnlen) :: cs,cs1,cs1colon,cs1plus,cs1times,cs2colon,cs2plus
162  character(len=fnlen) :: cs2times,cscolon,csplus,cstimes,trial_cs
163 !arrays
164  integer,allocatable :: int1(:),int2(:)
165  real(dp),allocatable :: dpr1(:),dpr2(:)
166 
167 ! *************************************************************************
168 
169  ABI_CHECK(marr >= narr, sjoin("marr", itoa(marr)," < narr ", itoa(narr)))
170 
171  ds_input_ = -1
172 
173  dozens=jdtset/10
174  unities=jdtset-10*dozens
175 
176  if(jdtset<0)then
177    write(message,'(a,i0,a)')' jdtset=',jdtset,', while it should be non-negative.'
178    MSG_ERROR(message)
179  end if
180 
181  if(jdtset>9999)then
182    write(message,'(a,i0,a)')' jdtset=',jdtset,', while it must be lower than 10000.'
183    MSG_ERROR(message)
184  end if
185 
186 !Default values : nothing has been read
187  itoken=0
188  opttoken=0
189 !Initialise flags in case of opttoken >= 2 later.
190  itoken_times=0
191  itoken_plus=0
192  itoken_colon=0
193  cslen=1
194 
195  if (narr/=0) then
196 
197    toklen=len_trim(token)
198 
199 !  --------------------------------------------------------------------------
200 !  (1) try to find the token with dataset number appended
201    if(jdtset>0)then
202 
203      call appdig(jdtset,'',appen)
204      cs=blank//token(1:toklen)//trim(appen)//blank
205      if(jdtset<10) then
206        cslen=toklen+3
207      else if(jdtset<100) then
208        cslen=toklen+4
209      else if(jdtset<1000) then
210        cslen=toklen+5
211      else if(jdtset<10000)then
212        cslen=toklen+6
213      end if
214 !    Map token to all upper case (make case-insensitive):
215      call inupper(cs)
216 !    Absolute index of blank//token//blank in string:
217      itoken=index(string,cs(1:cslen))
218 !    Look for another occurence of the same token in string, if so, leaves:
219      itoken2=index(string,cs(1:cslen), BACK=.true. )
220      if(itoken/=itoken2)then
221        write(message, '(7a)' )&
222 &       'There are two occurences of the keyword "',cs(1:cslen),'" in the input file.',ch10,&
223 &       'This is confusing, so it has been forbidden.',ch10,&
224 &       'Action: remove one of the two occurences.'
225        MSG_ERROR(message)
226      end if
227 
228      if(itoken/=0) then
229        opttoken=1
230        ds_input_=jdtset
231      end if
232    end if
233 
234 !  --------------------------------------------------------------------------
235 !  (2a) try to find the token appended with a string that contains the metacharacter "?".
236 
237    if(jdtset>0 .and. opttoken==0)then
238 
239 !    Use the metacharacter for the dozens, and save in cs and itoken
240      write(appen,'(i1)')unities
241      cs=blank//token(1:toklen)//'?'//trim(appen)//blank
242      cslen=toklen+4
243 !    Map token to all upper case (make case-insensitive):
244      call inupper(cs)
245 !    Absolute index of blank//token//blank in string:
246      itoken=index(string,cs(1:cslen))
247 !    Look for another occurence of the same token in string, if so, leaves:
248      itoken2=index(string,cs(1:cslen), BACK=.true. )
249      if(itoken/=itoken2)then
250        write(message, '(7a)' )&
251 &       'There are two occurences of the keyword "',cs(1:cslen),'" in the input file.',ch10,&
252 &       'This is confusing, so it has been forbidden.',ch10,&
253 &       'Action: remove one of the two occurences.'
254        MSG_ERROR(message)
255      end if
256      if(itoken/=0) then
257        opttoken=1
258        ds_input_=jdtset
259      end if
260 
261 !    Use the metacharacter for the unities, and save in cs1 and itoken1
262      write(appen,'(i1)')dozens
263      cs1=blank//token(1:toklen)//trim(appen)//'?'//blank
264 !    Map token to all upper case (make case-insensitive):
265      call inupper(cs1)
266 !    Absolute index of blank//token//blank in string:
267      itoken1=index(string,cs1(1:cslen))
268 !    Look for another occurence of the same token in string, if so, leaves:
269      itoken2=index(string,cs1(1:cslen), BACK=.true. )
270      if(itoken1/=itoken2)then
271        write(message, '(7a)' )&
272 &       'There are two occurences of the keyword "',cs1(1:cslen),'" in the input file.',ch10,&
273 &       'This is confusing, so it has been forbidden.',ch10,&
274 &       'Action: remove one of the two occurences.'
275        MSG_ERROR(message)
276      end if
277 
278      if(itoken/=0 .and. itoken1/=0)then
279        write(message, '(9a)' )&
280 &       'The keywords "',cs(1:cslen),'" and "',cs1(1:cslen),'"',ch10,&
281 &       'cannot be used together in the input file.',ch10,&
282 &       'Action: remove one of the two keywords.'
283        MSG_ERROR(message)
284      end if
285 
286      if(itoken1/=0)then
287        opttoken=1
288        itoken=itoken1
289        cs=cs1
290        ds_input_=jdtset
291      end if
292 
293    end if
294 
295 !  --------------------------------------------------------------------------
296 !  (2b) try to find the tokens defining a series
297    if(opttoken==0)then
298 
299      cs=token(1:toklen)
300 
301      cslen=toklen+3
302      cs1len=toklen+4
303 
304      cscolon=blank//token(1:toklen)//':'//blank
305      csplus=blank//token(1:toklen)//'+'//blank
306      cstimes=blank//token(1:toklen)//'*'//blank
307 
308      cs1colon=blank//token(1:toklen)//'?'//':'//blank
309      cs1plus=blank//token(1:toklen)//'?'//'+'//blank
310      cs1times=blank//token(1:toklen)//'?'//'*'//blank
311 
312      cs2colon=blank//token(1:toklen)//':'//'?'//blank
313      cs2plus=blank//token(1:toklen)//'+'//'?'//blank
314      cs2times=blank//token(1:toklen)//'*'//'?'//blank
315 
316 !    Map token to all upper case (make case-insensitive):
317      call inupper(cscolon)
318      call inupper(csplus)
319      call inupper(cstimes)
320      call inupper(cs1colon)
321      call inupper(cs1plus)
322      call inupper(cs1times)
323      call inupper(cs2colon)
324      call inupper(cs2plus)
325      call inupper(cs2times)
326 
327 !    Absolute index of tokens in string:
328      itoken_colon=index(string,cscolon(1:cslen))
329      itoken_plus=index(string,csplus(1:cslen))
330      itoken_times=index(string,cstimes(1:cslen))
331      itoken_1colon=index(string,cs1colon(1:cs1len))
332      itoken_1plus=index(string,cs1plus(1:cs1len))
333      itoken_1times=index(string,cs1times(1:cs1len))
334      itoken_2colon=index(string,cs2colon(1:cs1len))
335      itoken_2plus=index(string,cs2plus(1:cs1len))
336      itoken_2times=index(string,cs2times(1:cs1len))
337 
338 !    Look for another occurence of the same tokens in string
339      itoken2_colon=index(string,cscolon(1:cslen), BACK=.true. )
340      itoken2_plus=index(string,csplus(1:cslen), BACK=.true. )
341      itoken2_times=index(string,cstimes(1:cslen), BACK=.true. )
342      itoken2_1colon=index(string,cs1colon(1:cs1len), BACK=.true. )
343      itoken2_1plus=index(string,cs1plus(1:cs1len), BACK=.true. )
344      itoken2_1times=index(string,cs1times(1:cs1len), BACK=.true. )
345      itoken2_2colon=index(string,cs2colon(1:cs1len), BACK=.true. )
346      itoken2_2plus=index(string,cs2plus(1:cs1len), BACK=.true. )
347      itoken2_2times=index(string,cs2times(1:cs1len), BACK=.true. )
348 
349      if(jdtset==0)then
350 
351 !      If the multi-dataset mode is not used, no token should have been found
352        if(itoken_colon+itoken_plus+itoken_times+&
353 &       itoken_2colon+itoken_2plus+itoken_2times > 0 ) then
354          write(message,'(a,a,a,a,a,a,a,a,a,a,a,a, a)' )&
355 &         'Although the multi-dataset mode is not activated,',ch10,&
356 &         'the keyword "',trim(cs),'" has been found',ch10,&
357 &         'appended with  + * or :  .',ch10,&
358 &         'This is not allowed.',ch10,&
359 &         'Action: remove the appended keyword, or',ch10,&
360 &         'use the multi-dataset mode (ndtset/=0).'
361          MSG_ERROR(message)
362        end if
363        if(itoken_1colon+itoken_1plus+itoken_1times > 0 ) then
364          write(message, '(a,a,a,a,a,a,a,a,a,a,a,a,a)' )&
365 &         'Although the multi-dataset mode is not activated,',ch10,&
366 &         'the keyword "',trim(cs),'" has been found',ch10,&
367 &         'appended with ? , then + * or :  .',ch10,&
368 &         'This is not allowed.',ch10,&
369 &         'Action: remove the appended keyword, or',ch10,&
370 &         'use the multi-dataset mode (ndtset/=0).'
371          MSG_ERROR(message)
372        end if
373 
374      else
375 
376 !      If the multi-dataset mode is used, exactly zero or two token must be found
377        sum_token=0
378        if(itoken_colon/=0)sum_token=sum_token+1
379        if(itoken_plus /=0)sum_token=sum_token+1
380        if(itoken_times/=0)sum_token=sum_token+1
381        if(itoken_1colon/=0)sum_token=sum_token+1
382        if(itoken_1plus /=0)sum_token=sum_token+1
383        if(itoken_1times/=0)sum_token=sum_token+1
384        if(itoken_2colon/=0)sum_token=sum_token+1
385        if(itoken_2plus /=0)sum_token=sum_token+1
386        if(itoken_2times/=0)sum_token=sum_token+1
387 
388        if(sum_token/=0 .and. sum_token/=2) then
389          write(message, '(a,a,a,a,a,i3,a,a,a,a,a,a,a)' )&
390 &         'The keyword "',trim(cs),'" has been found to take part',ch10,&
391 &         'to series definition in the multi-dataset mode',sum_token,' times.',ch10,&
392 &         'This is not allowed, since it should be used once with ":",',ch10,&
393 &         'and once with "+" or "*".',ch10,&
394 &         'Action: change the number of occurences of this keyword.'
395          MSG_ERROR(message)
396        end if
397 
398 !      If the multi-dataset mode is used, make sure that
399 !      no twice the same combined keyword happens
400        ier=0
401        if(itoken_colon/=itoken2_colon)then
402          ier=1 ; cs=cscolon
403        end if
404        if(itoken_plus/=itoken2_plus)then
405          ier=1 ; cs=csplus
406        end if
407        if(itoken_times/=itoken2_times)then
408          ier=1 ; cs=cstimes
409        end if
410        if(itoken_1colon/=itoken2_1colon)then
411          ier=1 ; cs=cs1colon
412        end if
413        if(itoken_1plus/=itoken2_1plus)then
414          ier=1 ; cs=cs1plus
415        end if
416        if(itoken_1times/=itoken2_1times)then
417          ier=1 ; cs=cs1times
418        end if
419        if(itoken_2colon/=itoken2_2colon)then
420          ier=1 ; cs=cs2colon
421        end if
422        if(itoken_2plus/=itoken2_2plus)then
423          ier=1 ; cs=cs2plus
424        end if
425        if(itoken_2times/=itoken2_2times)then
426          ier=1 ; cs=cs2times
427        end if
428        if(ier==1)then
429          write(message, '(a,a,a,a,a,a,a)' )&
430 &         'There are two occurences of the keyword "',cs(1:cslen),'" in the input file.',ch10,&
431 &         'This is confusing, so it has been forbidden.',ch10,&
432 &         'Action: remove one of the two occurences.'
433          MSG_ERROR(message)
434        end if
435 
436 !      Select the series according to the presence of a colon flag
437        if(itoken_colon>0)then
438          opttoken=2
439          ds_input_=jdtset
440        else if(itoken_1colon>0)then
441          opttoken=3
442          cscolon=cs1colon ; csplus=cs1plus ; cstimes=cs1times
443          itoken_colon=itoken_1colon
444          itoken_plus=itoken_1plus ; itoken_times=itoken_1times
445          cslen=cs1len
446          ds_input_=jdtset
447        else if(itoken_2colon>0)then
448          opttoken=4
449          cscolon=cs2colon ; csplus=cs2plus ; cstimes=cs2times
450          itoken_colon=itoken_2colon
451          itoken_plus=itoken_2plus ; itoken_times=itoken_2times
452          cslen=cs1len
453          ds_input_=jdtset
454        end if
455 
456 !      Make sure that the proper combination of : + and * is found .
457        if(itoken_colon > 0 .and. (itoken_plus==0 .and. itoken_times==0) )then
458          write(message, '(13a)' )&
459 &         'The keyword "',cscolon(1:cslen),'" initiate a series,',ch10,&
460 &         'but there is no occurence of "',csplus(1:cslen),'" or "',cstimes(1:cslen),'".',ch10,&
461 &         'Action: either suppress the series, or make the increment',ch10,&
462 &         'or the factor available.'
463          MSG_ERROR(message)
464        end if
465        if(itoken_plus/=0 .and. itoken_times/=0)then
466          write(message, '(a,a, a,a,a,a,a)' )&
467 &         'The combined occurence of keywords "',csplus(1:cslen),'" and "',cstimes(1:cslen),'" is not allowed.',ch10,&
468 &         'Action: suppress one of them in your input file.'
469          MSG_ERROR(message)
470        end if
471        if(itoken_colon==0 .and. (itoken_plus/=0 .or. itoken_times/=0) ) then
472          cs=csplus
473          if(itoken_times/=0)cs=cstimes
474          write(message, '(a,a,a,a,a,a,a,a,a,a,a)' )&
475 &         'The keyword "',cscolon(1:cslen),'" does not appear in the input file.',ch10,&
476 &         'However, the keyword "',cs(1:cslen),'" appears.',ch10,&
477 &         'This is forbidden.',ch10,&
478 &         'Action: make the first appear, or suppress the second.'
479          MSG_ERROR(message)
480        end if
481 
482 !      At this stage, either
483 !      - itoken_colon vanish as well as itoken_plus and itoken_times
484 !      - itoken_colon does not vanish,
485 !      as well as one of itoken_plus or itoken_times
486 
487 !      End the condition of multi-dataset mode
488      end if
489 
490 !    End the check on existence of a series
491    end if
492 
493 !  --------------------------------------------------------------------------
494 !  (3) if not found, try to find the token with non-modified string
495    if(opttoken==0)then
496 
497      cs=blank//token(1:toklen)//blank
498      cslen=toklen+2
499 
500 !    Map token to all upper case (make case-insensitive):
501      call inupper(cs)
502 
503 !    Absolute index of blank//token//blank in string:
504      itoken=index(string,cs(1:cslen))
505 
506 !    Look for another occurence of the same token in string, if so, leaves:
507      itoken2=index(string,cs(1:cslen), BACK=.true. )
508      if(itoken/=itoken2)then
509        write(message, '(a,a,a,a,a,a,a)' )&
510 &       'There are two occurences of the keyword "',cs(1:cslen),'" in the input file.',ch10,&
511 &       'This is confusing, so it has been forbidden.',ch10,&
512 &       'Action: remove one of the two occurences.'
513        MSG_ERROR(message)
514      end if
515 
516      if(itoken/=0) then
517        opttoken=1
518        ds_input_=0
519      end if
520 
521    end if
522 
523 !  --------------------------------------------------------------------------
524 !  If jdtset==0, means that the multi-dataset mode is not used, so
525 !  checks whether the input file contains a multi-dataset keyword,
526 !  and if this occurs, stop. Check also the forbidden occurence of
527 !  use of 0 as a multi-dataset index.
528 !  Note that the occurence of series initiators has already been checked.
529 
530    do trial_jdtset=0,9
531      if(jdtset==0 .or. trial_jdtset==0)then
532        write(appen,'(i1)')trial_jdtset
533        trial_cs=blank//token(1:toklen)//trim(appen)
534        trial_cslen=toklen+2
535 !      Map token to all upper case (make case-insensitive):
536        call inupper(trial_cs)
537 !      Look for an occurence of this token in string, if so, leaves:
538        itoken2=index(string,trial_cs(1:trial_cslen))
539 !      If itoken2/=0
540        if(itoken2/=0)then
541          if(trial_jdtset==0)then
542            write(message, '(a,a,a,a,a,a,a)' )&
543 &           'There is an occurence of the keyword "',trim(token),'" appended with 0 in the input file.',ch10,&
544 &           'This is forbidden.',ch10,&
545 &           'Action: remove this occurence.'
546            call wrtout(std_out,message,'COLL')
547          else
548            write(message, '(a,a,a,a,a,i1,a,a,a,a,a)' )&
549 &           'In the input file, there is an occurence of the ',ch10,&
550 &           'keyword "',trim(token),'", appended with the digit "',trial_jdtset,'".',ch10,&
551 &           'This is forbidden when ndtset==0 .',ch10,&
552 &           'Action: remove this occurence, or change ndtset.'
553            call wrtout(std_out,message,'COLL')
554          end if
555          MSG_ERROR(message)
556        end if
557      end if
558    end do
559 
560  end if
561 
562 !===========================================================================
563 !At this stage, the location of the keyword string is known, as well
564 !as its length. So, can read the data.
565 !Usual reading if opttoken==1 (need itoken).
566 !If opttoken>=2, the characteristics of a series must be read
567 !(need itoken_colon and either itoken_plus or itoken_times)
568 
569  tread = 0
570  typevar='INT'
571  if(typevarphys=='LOG')typevar='INT'
572  if(typevarphys=='DPR' .or. typevarphys=='LEN' .or. typevarphys=='ENE' .or. typevarphys=='BFI')typevar='DPR'
573  if(typevarphys=='KEY')then
574    if(opttoken>=2)then
575      write(message, '(9a)' )&
576 &     'For the keyword "',cs(1:cslen),'", of KEY type,',ch10,&
577 &     'a series has been defined in the input file.',ch10,&
578 &     'This is forbidden.',ch10,&
579 &     'Action: check your input file.'
580      MSG_ERROR(message)
581    end if
582    if(narr>=2)then
583      write(message, '(9a)' )&
584 &     'For the keyword "',cs(1:cslen),'", of KEY type,',ch10,&
585 &     'the number of data requested is larger than 1.',ch10,&
586 &     'This is forbidden.',ch10,&
587 &     'Action: check your input file.'
588      MSG_ERROR(message)
589    end if
590    typevar='KEY'
591 !  write(std_out,*)' intagm : will read cs=',trim(cs)
592 !  stop
593  end if
594 
595 !There is something to be read if opttoken>=1
596  if(opttoken==1)then
597 
598 !  DEBUG
599 !  write(std_out,*)' intagm : opttoken==1 , token has been found, will read '
600 !  ENDDEBUG
601 
602 !  Absolute location in string of blank which follows token:
603    b1=itoken+cslen-1
604 
605 !  Read the array (or eventual scalar) that follows the blank
606 !  In case of typevarphys='KEY', the chain of character will be returned in cs.
607    call inarray(b1,cs,dprarr,intarr,marr,narr,string,typevarphys)
608 
609    if(typevarphys=='KEY')then
610      if (.not. PRESENT(key_value)) then
611        MSG_ERROR("typevarphys == KEY requires the optional argument key_value")
612      end if
613      !token=trim(cs)
614      !write(std_out,*)' intagm : after inarray, token=',trim(token)
615      key_value = TRIM(cs)
616    end if
617 
618 !  if this point is reached then data has been read in successfully
619    tread = 1
620 
621  else if(opttoken>=2)then
622 
623 !  write(std_out,*)' intagm : opttoken>=2 , token has been found, will read '
624    ABI_ALLOCATE(dpr1,(narr))
625    ABI_ALLOCATE(dpr2,(narr))
626    ABI_ALLOCATE(int1,(narr))
627    ABI_ALLOCATE(int2,(narr))
628 
629 !  Absolute location in string of blank which follows token//':':
630    b1=itoken_colon+cslen-1
631    call inarray(b1,cscolon,dpr1,int1,narr,narr,string,typevarphys)
632 
633 !  Initialise number even if the if series treat all cases.
634    number=1
635 !  Define the number of the term in the series
636    if(opttoken==2)number=jdtset-1
637    if(opttoken==3)number=unities-1
638    if(opttoken==4)number=dozens-1
639 
640 !  Distinguish additive and multiplicative series
641    if(itoken_plus/=0)then
642 
643      b1=itoken_plus+cslen-1
644      call inarray(b1,csplus,dpr2,int2,narr,narr,string,typevarphys)
645 
646      if(typevar=='INT')then
647        intarr(1:narr)=int1(:)+int2(:)*number
648      else if(typevar=='DPR')then
649        dprarr(1:narr)=dpr1(:)+dpr2(:)*number
650      end if
651 
652    else if(itoken_times/=0)then
653 
654      b1=itoken_times+cslen-1
655      call inarray(b1,cstimes,dpr2,int2,narr,narr,string,typevarphys)
656      if(typevar=='INT')then
657        intarr(1:narr)=int1(:)*int2(:)**number
658      else if(typevar=='DPR')then
659        dprarr(1:narr)=dpr1(:)*dpr2(:)**number
660      end if
661 
662    end if
663 
664    tread = 1
665 
666    ABI_DEALLOCATE(dpr1)
667    ABI_DEALLOCATE(dpr2)
668    ABI_DEALLOCATE(int1)
669    ABI_DEALLOCATE(int2)
670  end if
671 
672  if(present(ds_input)) then
673    ds_input = ds_input_
674  end if
675 
676 !DEBUG
677 !write(std_out,*) ' intagm : exit value tread=',tread
678 !write(std_out,*) ' intarr =',intarr(1:narr)
679 !write(std_out,*) ' dprarr =',dprarr(1:narr)
680 !stop
681 !ENDDEBUG
682 
683 end subroutine intagm