TABLE OF CONTENTS


ABINIT/ingeobld [ Functions ]

[ Top ] [ Functions ]

NAME

 ingeobld

FUNCTION

 The geometry builder.
 Start from the types and coordinates of the primitive atoms
 and produce the completed set of atoms, by using the definition
 of objects, then application of rotation, translation and repetition.

COPYRIGHT

 Copyright (C) 1998-2018 ABINIT group (XG)
 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

 iout=unit number of output file
 jdtset=number of the dataset looked for
 lenstr=actual length of the string
 natrd=number of atoms that have been read in the calling routine
 natom=number of atoms
 nobj=the number of objects
 string*(*)=character string containing all the input data, used
  only if choice=1 or 3. Initialized previously in instrng.
 typat_read(natrd)=type integer for each atom in the primitive set
 xcart_read(3,natrd)=cartesian coordinates of atoms (bohr), in the primitive set

OUTPUT

 typat(natom)=type integer for each atom in cell
 xcart(3,natom)=cartesian coordinates of atoms (bohr)

PARENTS

      ingeo

CHILDREN

      intagm,wrtout

SOURCE

 44 #if defined HAVE_CONFIG_H
 45 #include "config.h"
 46 #endif
 47 
 48 #include "abi_common.h"
 49 
 50 subroutine ingeobld (iout,jdtset,lenstr,natrd,natom,nobj,string,typat,typat_read,xcart,xcart_read)
 51 
 52  use defs_basis
 53  use m_profiling_abi
 54  use m_errors
 55 
 56  use m_parser,  only : intagm
 57 
 58 !This section has been created automatically by the script Abilint (TD).
 59 !Do not modify the following lines by hand.
 60 #undef ABI_FUNC
 61 #define ABI_FUNC 'ingeobld'
 62  use interfaces_14_hidewrite
 63 !End of the abilint section
 64 
 65  implicit none
 66 
 67 !Arguments ------------------------------------
 68 !scalars
 69  integer,intent(in) :: iout,jdtset,lenstr,natom,natrd,nobj
 70  character(len=*),intent(in) :: string
 71 !arrays
 72  integer,intent(in) :: typat_read(natrd)
 73  integer,intent(out) :: typat(natom)
 74  real(dp),intent(in) :: xcart_read(3,natrd)
 75  real(dp),intent(out) :: xcart(3,natom)
 76 
 77 !Local variables-------------------------------
 78  character(len=*), parameter :: format01110 ="(1x,a6,1x,(t9,8i8) )"
 79  character(len=*), parameter :: format01160 ="(1x,a6,1x,1p,(t9,3g18.10)) "
 80 !scalars
 81  integer :: belonga,belongb,iatom,iatrd,ii,irep,irep1,irep2,irep3,ivac,marr
 82  integer :: natom_toberead,nread,objan,objbn,rotate,shift,tread,vacnum
 83  real(dp) :: angle,cosine,norm2per,norma,normb,normper,project,sine
 84  character(len=500) :: message
 85 !arrays
 86  integer :: objarf(3),objbrf(3)
 87  integer,allocatable :: objaat(:),objbat(:),vaclst(:)
 88  real(dp) :: axis2(3),axis3(3),axisa(3),axisb(3),objaax(6),objaro(4),objatr(12)
 89  real(dp) :: objbax(6),objbro(4),objbtr(12),parall(3),perpen(3),rotated(3)
 90  real(dp) :: vectora(3),vectorb(3)
 91  real(dp),allocatable :: typat_full(:),xcart_full(:,:)
 92 !no_abirules
 93 !Dummy arguments for subroutine 'intagm' to parse input file
 94  integer,allocatable :: intarr(:)
 95  real(dp),allocatable :: dprarr(:)
 96 
 97 ! *************************************************************************
 98 
 99  marr=max(12,3*natom)
100  ABI_ALLOCATE(intarr,(marr))
101  ABI_ALLOCATE(dprarr,(marr))
102 
103 !1) Set up the number of vacancies.
104 
105 !This is the default
106  vacnum=0
107  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'vacnum',tread,'INT')
108  if(tread==1) vacnum=intarr(1)
109 
110  if (vacnum>0)then
111    ABI_ALLOCATE(vaclst,(vacnum))
112 !  Read list of atoms to be suppressed to create vacancies
113    call intagm(dprarr,intarr,jdtset,marr,vacnum,string(1:lenstr),'vaclst',tread,'INT')
114    if(tread==1) vaclst(:)=intarr(1:vacnum)
115    if(tread/=1)then
116      write(message, '(a,a,a,a,a)' )&
117 &     'The array vaclst MUST be initialized in the input file',ch10,&
118 &     'when vacnum is non-zero.',ch10,&
119 &     'Action: initialize vaclst in your input file.'
120      MSG_ERROR(message)
121    end if
122  end if
123 
124  natom_toberead=natom+vacnum
125 
126 !2) Set up list and number of atoms in objects, and the --------------
127 !operations to be performed on objects.
128 
129  write(message,'(80a,a)')('=',ii=1,80),ch10
130  call wrtout(std_out,message,'COLL')
131  call wrtout(iout,message,'COLL')
132 
133  write(message, '(a,a)' )&
134 & '--ingeobld: echo values of variables connected to objects --------',ch10
135  call wrtout(std_out,message,'COLL')
136  call wrtout(iout,message,'COLL')
137 
138  if(vacnum>0)then
139    write(iout,format01110) 'vacnum',vacnum
140    write(std_out,format01110) 'vacnum',vacnum
141    write(iout,'(1x,a6,1x,(t9,20i3))') 'vaclst',vaclst(:)
142    write(std_out,'(1x,a6,1x,(t9,20i3))') 'vaclst',vaclst(:)
143    write(iout, '(a)' ) ' '
144    write(std_out,'(a)' ) ' '
145  end if
146 
147  write(iout,format01110) 'nobj',nobj
148  write(std_out,format01110) 'nobj',nobj
149 
150  if(nobj/=1 .and. nobj/=2)then
151    write(message, '(a,a,a,i8,a,a,a)' )&
152 &   'The number of object (nobj) must be either 1 or 2,',ch10,&
153 &   'while the input file has  nobj=',nobj,'.',ch10,&
154 &   'Action: correct nobj in your input file.'
155    MSG_ERROR(message)
156  end if
157 
158  if(nobj==1 .or. nobj==2)then
159 
160 !  Read the number of atoms of the object a
161    call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'objan',tread,'INT')
162    if(tread==1) objan=intarr(1)
163 
164    if(tread/=1)then
165      write(message, '(a,a,a,i8,a,a,a,a,a)' )&
166 &     'The number of atoms in object a (objan) must be initialized',ch10,&
167 &     'in the input file, when nobj=',nobj,'.',ch10,&
168 &     'This is not the case.',ch10,&
169 &     'Action: correct objan in your input file.'
170      MSG_ERROR(message)
171    end if
172 
173    write(iout, '(a)' ) ' '
174    write(std_out,'(a)' ) ' '
175    write(iout,format01110) 'objan',objan
176    write(std_out,format01110) 'objan',objan
177 
178    if(objan<=1 .or. objan>natom)then
179      write(message, '(a,a,a,a,a,i8,a,a,a)' )&
180 &     'The number of atoms in object a (objan) must be larger than 0',ch10,&
181 &     'and smaller than natom.',ch10,&
182 &     'It is equal to ',objan,', an unacceptable value.',ch10,&
183 &     'Action: correct objan in your input file.'
184      MSG_ERROR(message)
185    end if
186 
187 !  Read list of atoms in object a
188    call intagm(dprarr,intarr,jdtset,marr,objan,string(1:lenstr),'objaat',tread,'INT')
189    ABI_ALLOCATE(objaat,(objan))
190    if(tread==1) objaat(1:objan)=intarr(1:objan)
191 
192    if(tread/=1)then
193      write(message, '(a,a,a,i8,a,a,a,a,a)' )&
194 &     'The list of atoms in object a (objaat) must be initialized',ch10,&
195 &     'in the input file, when nobj=',nobj,'.',ch10,&
196 &     'This is not the case.',ch10,&
197 &     'Action: initialize objaat in your input file.'
198      MSG_ERROR(message)
199    end if
200 
201    write(iout,'(1x,a6,1x,(t9,20i3))') 'objaat',objaat(:)
202    write(std_out,'(1x,a6,1x,(t9,20i3))') 'objaat',objaat(:)
203 
204    do iatom=1,objan
205      if(objaat(iatom)<1 .or. objaat(iatom)>natom)then
206        write(message, '(a,i8,a,a,i8,4a)' )&
207 &       'The input value of objaat for atom number ',iatom,ch10,&
208 &       'is equal to ',objaat(iatom),', an unacceptable value :',ch10,&
209 &       'it should be between 1 and natom. ',&
210 &       'Action: correct the array objaat in your input file.'
211        MSG_ERROR(message)
212      end if
213    end do
214 
215    if(objan>1)then
216      do iatom=1,objan-1
217        if( objaat(iatom)>=objaat(iatom+1) )then
218          write(message, '(a,i8,a,a,a,a,a,a)' )&
219 &         'The input value of objaat for atom number ',iatom,ch10,&
220 &         'is larger or equal to the one of the next atom,',ch10,&
221 &         'while this list should be ordered, and an atom cannot be repeated.',ch10,&
222 &         'Action: correct the array objaat in your input file.'
223          MSG_ERROR(message)
224        end if
225      end do
226    end if
227 
228 !  Read repetition factors
229    objarf(1:3)=1
230    call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'objarf',tread,'INT')
231    if(tread==1) objarf(1:3)=intarr(1:3)
232    write(iout,'(1x,a6,1x,(t9,20i3))') 'objarf',objarf(:)
233    write(std_out,'(1x,a6,1x,(t9,20i3))') 'objarf',objarf(:)
234 
235    if(tread==1)then
236      do irep=1,3
237        if(objarf(irep)<1)then
238          write(message, '(a,a,a,3i8,a,a,a)' )&
239 &         'The input values of objarf(1:3) must be positive,',ch10,&
240 &         'while it is ',objarf(1:3),'.',ch10,&
241 &         'Action: correct objarf in your input file.'
242          MSG_ERROR(message)
243        end if
244      end do
245    end if
246 
247 !  Modify the number of atoms to be read
248    natom_toberead=natom_toberead-objan*(objarf(1)*objarf(2)*objarf(3)-1)
249 
250 !  Read rotations angles and translations
251    objaro(1:4)=0.0_dp
252    objatr(1:12)=0.0_dp
253    if (objarf(1)*objarf(2)*objarf(3) ==1) then
254      nread=1
255    else if (objarf(2)*objarf(3) ==1) then
256      nread=2
257    else if (objarf(3) ==1) then
258      nread=3
259    else
260      nread=4
261    end if
262    call intagm(dprarr,intarr,jdtset,marr,nread,string(1:lenstr),'objaro',tread,'DPR')
263    if(tread==1) objaro(1:nread)=dprarr(1:nread)
264 
265    call intagm(dprarr,intarr,jdtset,marr,3*nread,string(1:lenstr),'objatr',tread,'LEN')
266 
267    if(tread==1) objatr(1:3*nread)=dprarr(1:3*nread)
268    write(iout,format01160) 'objaro',objaro(1:4)
269    write(std_out,format01160) 'objaro',objaro(1:4)
270    write(iout,format01160) 'objatr',objatr(1:12)
271    write(std_out,format01160) 'objatr',objatr(1:12)
272 !  If needed, read axes, but default to the x-axis to avoid errors later
273    objaax(1:6)=0.0_dp ; objaax(4)=1.0_dp
274 
275    if(abs(objaro(1))+abs(objaro(2))+abs(objaro(3))+abs(objaro(4)) > 1.0d-10) then
276      call intagm(dprarr,intarr,jdtset,marr,6,string(1:lenstr),'objaax',tread,'LEN')
277      if(tread==1) objaax(1:6)=dprarr(1:6)
278      if(tread/=1)then
279        write(message, '(a,a,a,a,a,a,a)' )&
280 &       'The axis of object a (objaax) must be initialized',ch10,&
281 &       'in the input file, when rotations (objaro) are present.',ch10,&
282 &       'This is not the case.',ch10,&
283 &       'Action: initialize objaax in your input file.'
284        MSG_ERROR(message)
285      end if
286      write(iout,format01160) 'objaax',objaax(1:6)
287      write(std_out,format01160) 'objaax',objaax(1:6)
288    end if
289 
290    axisa(1:3)=objaax(4:6)-objaax(1:3)
291    norma=axisa(1)**2+axisa(2)**2+axisa(3)**2
292 
293    if(norma<1.0d-10)then
294      write(message, '(5a)' )&
295 &     'The two points defined by the input array objaax are too',ch10,&
296 &     'close to each other, and will not be used to define an axis.',ch10,&
297 &     'Action: correct objaax in your input file.'
298      MSG_ERROR(message)
299    end if
300    axisa(1:3)=axisa(1:3)/sqrt(norma)
301 
302 !  End condition of existence of a first object
303  end if
304 
305  if(nobj==2)then
306 
307 !  Read the number of atoms of the object b
308    call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'objbn',tread,'INT')
309    if(tread==1) objbn=intarr(1)
310 
311    if(tread/=1)then
312      write(message, '(a,a,a,i8,a,a,a,a,a)' )&
313 &     'The number of atoms in object b (objbn) must be initialized',ch10,&
314 &     'in the input file, when nobj=',nobj,'.',ch10,&
315 &     'This is not the case.',ch10,&
316 &     'Action: initialize objbn in your input file.'
317      MSG_ERROR(message)
318    end if
319 
320    write(iout, '(a)' ) ' '
321    write(std_out,'(a)' ) ' '
322    write(iout,format01110) 'objbn',objbn
323    write(std_out,format01110) 'objbn',objbn
324 
325    if(objbn<=1 .or. objbn>natom)then
326      write(message, '(a,a,a,a,a,i8,a,a,a)' )&
327 &     'The number of atoms in object b (objbn) must be larger than 0',ch10,&
328 &     'and smaller than natom.',ch10,&
329 &     'It is equal to ',objbn,', an unacceptable value.',ch10,&
330 &     'Action: correct objbn in your input file.'
331      MSG_ERROR(message)
332    end if
333 
334 !  Read list of atoms in object b
335    call intagm(dprarr,intarr,jdtset,marr,objbn,string(1:lenstr),'objbat',tread,'INT')
336    ABI_ALLOCATE(objbat,(objbn))
337 
338    if(tread==1) objbat(1:objbn)=intarr(1:objbn)
339    if(tread/=1)then
340      write(message, '(a,a,a,i8,a,a,a,a,a)' )&
341 &     'The list of atoms in object b (objbat) must be initialized',ch10,&
342 &     'in the input file, when nobj=',nobj,'.',ch10,&
343 &     'This is not the case.',ch10,&
344 &     'Action: initialize objbat in your input file.'
345      MSG_ERROR(message)
346    end if
347 
348    write(iout,'(1x,a6,1x,(t9,20i3))') 'objbat',objbat(:)
349    write(std_out,'(1x,a6,1x,(t9,20i3))') 'objbat',objbat(:)
350 
351    do iatom=1,objbn
352      if(objbat(iatom)<1 .or. objbat(iatom)>natom)then
353        write(message, '(a,i8,a,a,i8,a,a,a,a,a)' )&
354 &       'The input value of objbat for atom number ',iatom,ch10,&
355 &       'is equal to ',objbat(iatom),', an unacceptable value :',ch10,&
356 &       'it should be between 1 and natom. ',ch10,&
357 &       'Action: correct objbat in your input file.'
358        MSG_ERROR(message)
359      end if
360    end do
361 
362    if(objbn>1)then
363      do iatom=1,objbn-1
364        if( objbat(iatom)>=objbat(iatom+1) )then
365          write(message, '(a,i8,a,a,a,a,a,a)' )&
366 &         'The input value of objbat for atom number ',iatom,ch10,&
367 &         'is larger or equal to the one of the next atom,',ch10,&
368 &         'while this list should be ordered, and an atom cannot be repeated.',ch10,&
369 &         'Action: correct the array objbat in the input file.'
370          MSG_ERROR(message)
371        end if
372      end do
373    end if
374 
375 !  Read repetition factors
376    objbrf(1:3)=1
377    call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'objbrf',tread,'INT')
378    if(tread==1) objbrf(1:3)=intarr(1:3)
379    write(iout,'(1x,a6,1x,(t9,20i3))') 'objbrf',objbrf(:)
380    write(std_out,'(1x,a6,1x,(t9,20i3))') 'objbrf',objbrf(:)
381 
382    if(tread==1)then
383      do irep=1,3
384        if(objbrf(irep)<1)then
385          write(message, '(a,a,a,3i8,a,a,a)' )&
386 &         'The input values of objbrf(1:3) must be positive,',ch10,&
387 &         'while it is ',objbrf(1:3),'.',ch10,&
388 &         'Action: correct objbrf in your input file.'
389          MSG_ERROR(message)
390        end if
391      end do
392    end if
393 
394 !  Modify the number of atoms to be read
395    natom_toberead=natom_toberead-objbn*(objbrf(1)*objbrf(2)*objbrf(3)-1)
396 !  Read rotations angles and translations
397    objbro(1:4)=0.0_dp
398    objbtr(1:12)=0.0_dp
399    if (objbrf(1)*objbrf(2)*objbrf(3) ==1) then
400      nread=1
401    else if (objbrf(2)*objbrf(3) ==1) then
402      nread=2
403    else if (objbrf(3) ==1) then
404      nread=3
405    else
406      nread=4
407    end if
408    call intagm(dprarr,intarr,jdtset,marr,nread,string(1:lenstr),'objbro',tread,'DPR')
409    if(tread==1) objbro(1:nread)=dprarr(1:nread)
410 
411    call intagm(dprarr,intarr,jdtset,marr,3*nread,string(1:lenstr),'objbtr',tread,'LEN')
412    if(tread==1) objbtr(1:3*nread)=dprarr(1:3*nread)
413 
414    write(iout,format01160) 'objbro',objbro(1:4)
415    write(std_out,format01160) 'objbro',objbro(1:4)
416    write(iout,format01160) 'objbtr',objbtr(1:12)
417    write(std_out,format01160) 'objbtr',objbtr(1:12)
418 
419 !  If needed, read axes, but default to the x-axis to avoid errors later
420    objbax(1:6)=0.0_dp ; objbax(4)=1.0_dp
421    if(abs(objbro(1))+abs(objbro(2))+abs(objbro(3))+abs(objbro(4)) > 1.0d-10) then
422      call intagm(dprarr,intarr,jdtset,marr,6,string(1:lenstr),'objbax',tread,'LEN')
423      if(tread==1) objbax(1:6)=dprarr(1:6)
424      if(tread/=1)then
425        write(message, '(a,a,a,a,a,a,a)' )&
426 &       'The axis of object b (objbax) must be initialized',ch10,&
427 &       'in the input file, when rotations (objbro) are present.',ch10,&
428 &       'This is not the case.',ch10,&
429 &       'Action: initialize objbax in your input file.'
430        MSG_ERROR(message)
431      end if
432      write(iout,format01160) 'objbax',objbax(1:6)
433      write(std_out,format01160) 'objbax',objbax(1:6)
434    end if
435    axisb(1:3)=objbax(4:6)-objbax(1:3)
436    normb=axisb(1)**2+axisb(2)**2+axisb(3)**2
437    if(normb<1.0d-10)then
438      write(message, '(5a)' )&
439 &     'The two points defined by the input array objbax are too',ch10,&
440 &     'close to each other, and will not be used to define an axis.',ch10,&
441 &     'Action: correct objbax in your input file.'
442      MSG_ERROR(message)
443    end if
444    axisb(1:3)=axisb(1:3)/sqrt(normb)
445 
446 !  Check whether both lists are disjoints. Use a very primitive algorithm.
447    do iatom=1,objan
448      do ii=1,objbn
449        if(objaat(iatom)==objbat(ii))then
450          write(message, '(6a,i8,a,i8,3a)' )&
451 &         'The objects a and b cannot have a common atom, but it is',ch10,&
452 &         'found that the values of objaat and objbat ',&
453 &         ' are identical, for their',ch10,&
454 &         'atoms number ',iatom,' and ',ii,'.',ch10,&
455 &         'Action: change objaat and/or objbat so that they have no common atom anymore.'
456          MSG_ERROR(message)
457        end if
458      end do
459    end do
460 
461 !  End condition of existence of a second object
462  end if
463 
464 !Check whether the number of atoms to be read obtained by relying
465 !on natom, vacnum and the object definitions, or from natrd coincide
466  if(natrd/=natom_toberead)then
467    write(message,'(11a,i0,a,i0,2a,i0,a)' )&
468 &   ' ingeobld : ERROR -- ',ch10,&
469 &   '  The number of atoms to be read (natrd) must be equal',ch10,&
470 &   '  to the total number of atoms (natom), plus',ch10,&
471 &   '  the number of vacancies (vacnum), minus',ch10,&
472 &   '  the number of atoms added by the repetition of objects.',ch10,&
473 &   '  This is not the case : natrd= ',natrd,', natom= ',natom,ch10,&
474 &   ', vacnum= ',vacnum,';'
475    call wrtout(std_out,message,"COLL")
476 
477    if(nobj==1 .or. nobj==2) then
478      write(message,'(a,i3,a,3i3,a,i5,a)' )&
479 &     '   object a : objan=',objan,', objarf(1:3)=',objarf(1:3),&
480 &     ' => adds ',objan*(objarf(1)*objarf(2)*objarf(3)-1),' atoms.'
481      call wrtout(std_out,message,"COLL")
482    end if
483 
484    if(nobj==2) then
485      write(message,'(a,i3,a,3i3,a,i5,a)' )&
486 &     '   object b : objbn=',objbn,', objbrf(1:3)=',objbrf(1:3),&
487 &     ' => adds ',objbn*(objbrf(1)*objbrf(2)*objbrf(3)-1),' atoms.'
488      call wrtout(std_out,message,"COLL")
489    end if
490 
491    write(message,'(3a)' )&
492 &   '  Action : check the correspondence between natom+vacnum on one side,',ch10,&
493 &   '           and natrd, objan, objbn, objarf and objbrf on the other side.'
494    MSG_ERROR(message)
495  end if
496 
497 !6) Produce full set of atoms
498 
499 !Print the initial atom coordinates if the geometry builder is used
500  write(iout, '(/,a)' )  ' Cartesian coordinates of the primitive atoms '
501  write(std_out,'(/,a)' )' Cartesian coordinates of the primitive atoms '
502  write(iout,format01160) '      ',xcart_read(:,:)
503  write(std_out,format01160) '      ',xcart_read(:,:)
504 
505  ABI_ALLOCATE(typat_full,(natom+vacnum))
506  ABI_ALLOCATE(xcart_full,(3,natom+vacnum))
507 
508 !Use the work array xcart_full to produce full set of atoms,
509 !including those coming from repeated objects.
510  iatom=1
511  do iatrd=1,natrd
512 
513    belonga=0 ; belongb=0
514    if(nobj==1 .or. nobj==2)then
515 !    Determine whether the atom belongs to object a
516      do ii=1,objan
517        if(iatrd==objaat(ii))belonga=ii
518      end do
519    end if
520    if(nobj==2)then
521 !    Determine whether the atom belong to object b
522      do ii=1,objbn
523        if(iatrd==objbat(ii))belongb=ii
524      end do
525    end if
526 
527    write(std_out,'(a,i5,a,i2,i2,a)' ) &
528 &   ' ingeobld : treating iatrd=',iatrd,', belong(a,b)=',belonga,belongb,'.'
529 
530 !  In case it does not belong to an object
531    if(belonga==0 .and. belongb==0)then
532      xcart_full(1:3,iatom)=xcart_read(1:3,iatrd)
533      typat_full(iatom)=typat_read(iatrd)
534      iatom=iatom+1
535    else
536 
537 !    Repeat, rotate and translate this atom
538      if(belonga/=0)then
539 
540 !      Treat object a
541 !      Compute the relative coordinate of atom with respect to first point of axis
542        vectora(1:3)=xcart_read(1:3,iatrd)-objaax(1:3)
543 !      Project on axis
544        project=vectora(1)*axisa(1)+vectora(2)*axisa(2)+vectora(3)*axisa(3)
545 !      Get the parallel part
546        parall(1:3)=project*axisa(1:3)
547 !      Get the perpendicular part, to be rotated
548        perpen(1:3)=vectora(1:3)-parall(1:3)
549 !      Compute the norm of the perpendicular part
550        norm2per=perpen(1)**2+perpen(2)**2+perpen(3)**2
551 !      Initialisation to avoid warnings even if used behind if rotate == 1.
552        normper = 0
553 !      It the norm is too small, there is not need to rotate
554        rotate=0
555        if(norm2per>=1.0d-18)then
556          rotate=1
557          normper=sqrt(norm2per)
558          axis2(1:3)=perpen(1:3)/normper
559 !        Get the vector perpendicular to axisa and axisa2
560          axis3(1)=axisa(2)*axis2(3)-axisa(3)*axis2(2)
561          axis3(2)=axisa(3)*axis2(1)-axisa(1)*axis2(3)
562          axis3(3)=axisa(1)*axis2(2)-axisa(2)*axis2(1)
563        end if
564 
565 !      Here the repetition loop
566        do irep3=1,objarf(3)
567          do irep2=1,objarf(2)
568            do irep1=1,objarf(1)
569 !            Here the rotation
570              if(rotate==1)then
571 !              Compute the angle of rotation
572                angle=objaro(1)+(irep1-1)*objaro(2)+                     &
573 &               (irep2-1)*objaro(3)+(irep3-1)*objaro(4)
574                cosine=cos(angle/180.0*pi)
575                sine=sin(angle/180.0*pi)
576                rotated(1:3)=objaax(1:3)+parall(1:3)+&
577 &               normper*(cosine*axis2(1:3)+sine*axis3(1:3))
578              else
579                rotated(1:3)=vectora(1:3)
580              end if
581 !            Here the translation
582              xcart_full(1:3,iatom)=rotated(1:3)+objatr(1:3)+&
583 &             (irep1-1)*objatr(4:6)+(irep2-1)*objatr(7:9)+(irep3-1)*objatr(10:12)
584              typat_full(iatom)=typat_read(iatrd)
585              iatom=iatom+1
586            end do
587          end do
588 !        End the repetition loop
589        end do
590 
591      else
592 !      If the atom belong to object b
593 !      Compute the relative coordinate of atom with respect to first point of axis
594        vectorb(1:3)=xcart_read(1:3,iatrd)-objbax(1:3)
595 !      Project on axis
596        project=vectorb(1)*axisb(1)+vectorb(2)*axisb(2)+vectorb(3)*axisb(3)
597 !      Get the parallel part
598        parall(1:3)=project*axisb(1:3)
599 !      Get the perpendicular part, to be rotated
600        perpen(1:3)=vectorb(1:3)-parall(1:3)
601 !      Compute the norm of the perpendicular part
602        norm2per=perpen(1)**2+perpen(2)**2+perpen(3)**2
603 !      Initialisation to avoid warnings even if used behind if rotate == 1.
604        normper = 0
605 !      It the norm is too small, there is not need to rotate
606        rotate=0
607        if(norm2per>=1.0d-18)then
608          rotate=1
609          normper=sqrt(norm2per)
610          axis2(1:3)=perpen(1:3)/normper
611 !        Get the vector perpendicular to axisb and axis2
612          axis3(1)=axisb(2)*axis2(3)-axisb(3)*axis2(2)
613          axis3(2)=axisb(3)*axis2(1)-axisb(1)*axis2(3)
614          axis3(3)=axisb(1)*axis2(2)-axisb(2)*axis2(1)
615        end if
616 !      Here the repetition loop
617        do irep3=1,objbrf(3)
618          do irep2=1,objbrf(2)
619            do irep1=1,objbrf(1)
620 !            Here the rotation
621              if(rotate==1)then
622 !              Compute the angle of rotation
623                angle=objbro(1)+(irep1-1)*objbro(2)+                      &
624 &               (irep2-1)*objbro(3)+ (irep3-1)*objbro(4)
625                cosine=cos(angle/180.0*pi)
626                sine=sin(angle/180.0*pi)
627                rotated(1:3)=objbax(1:3)+parall(1:3)+&
628 &               normper*(cosine*axis2(1:3)+sine*axis3(1:3))
629              else
630                rotated(1:3)=vectorb(1:3)
631              end if
632 !            Here the translation
633              xcart_full(1:3,iatom)=rotated(1:3)+objbtr(1:3)+&
634 &             (irep1-1)*objbtr(4:6)+(irep2-1)*objbtr(7:9)+(irep3-1)*objbtr(10:12)
635              typat_full(iatom)=typat_read(iatrd)
636              iatom=iatom+1
637            end do
638          end do
639 !        End the repetition loop
640        end do
641 
642 !      End the condition of belonging to object b
643      end if
644 
645 !    End the condition of belonging to an object
646    end if
647 
648 !  End the loop on atoms
649  end do
650 
651 !Create the vacancies here
652  if(vacnum/=0)then
653 !  First label the vacant atoms as belonging to typat 0
654    do ivac=1,vacnum
655      typat_full(vaclst(ivac))=0
656    end do
657 !  Then compact the arrays
658    shift=0
659    do iatom=1,natom
660      if(typat_full(iatom+shift)==0) shift=shift+1
661      if(shift/=0)then
662        xcart_full(1:3,iatom)=xcart_full(1:3,iatom+shift)
663        typat_full(iatom)=typat_full(iatom+shift)
664      end if
665    end do
666  end if
667 
668 !Transfer the content of xcart_full and typat_full to the proper
669 !location
670  xcart(:,1:natom)=xcart_full(:,1:natom)
671  typat(1:natom)=typat_full(1:natom)
672 
673  ABI_DEALLOCATE(typat_full)
674  ABI_DEALLOCATE(xcart_full)
675  if(allocated(objaat)) then
676    ABI_DEALLOCATE(objaat)
677  end if
678  if(allocated(objbat)) then
679    ABI_DEALLOCATE(objbat)
680  end if
681 
682  ABI_DEALLOCATE(intarr)
683  ABI_DEALLOCATE(dprarr)
684  if (vacnum>0)  then
685    ABI_DEALLOCATE(vaclst)
686  end if
687 
688 end subroutine ingeobld