1*
2* $Id$
3*
4
5*     ***************************
6*     *             		*
7*     *        cpsi_write	*
8*     *             		*
9*     ***************************
10
11      subroutine cpsi_write(spin_orbit,ispin,ne,nbrillioun,psi2_tag)
12      implicit none
13      logical spin_orbit
14      integer ispin,ne(2),nbrillioun
15      integer psi2_tag
16
17#include "bafdecls.fh"
18#include "errquit.fh"
19#include "stdio.fh"
20#include "util.fh"
21
22*    *** local variables ***
23      integer version,l,ispin0
24      integer nfft3d
25      integer nfft(3)
26      real*8  unita(3,3)
27      character*50 filename
28      character*255 full_filename
29
30      integer MASTER,taskid,taskid_j,taskid_k
31      parameter(MASTER=0)
32      integer n,nb,indx,occupation
33      integer q,pj,nbq,pk
34
35c     complex*16 tmp(*)
36      integer tmp(2),tmp2(2),tmpe(2)
37      logical value
38
39*     ***** local functions ****
40      logical      control_print
41      external     control_print
42      character*50 control_output_psi
43      external     control_output_psi
44      double precision control_unita
45      external         control_unita
46      integer  control_ngrid,cpsi_data_get_chnk,cpsi_data_get_next
47      integer  cpsi_data_get_ptr
48      external control_ngrid,cpsi_data_get_chnk,cpsi_data_get_next
49      external cpsi_data_get_ptr
50
51      call Parallel_taskid(taskid)
52      call Parallel3d_taskid_j(taskid_j)
53      call Parallel3d_taskid_k(taskid_k)
54      call C3dB_nfft3d(1,nfft3d)
55      occupation = cpsi_data_get_next(psi2_tag)
56
57      value = BA_push_get(mt_dcpl,nfft3d,'tmp',tmp(2),tmp(1))
58      if (.not. value) call errquit('out of stack memory',0, MA_ERR)
59      value = BA_push_get(mt_dcpl,nfft3d,'tmp2',tmp2(2),tmp2(1))
60      if (.not. value) call errquit('out of stack memory',0, MA_ERR)
61      value = BA_push_get(mt_dbl,(ne(1)+ne(2)),'tmpe',tmpe(2),tmpe(1))
62      if (.not. value) call errquit('out of stack memory',0, MA_ERR)
63
64      version = 5 ! version 5 - band, version 3 and 4 - pspw
65      nfft(1) = control_ngrid(1)
66      nfft(2) = control_ngrid(2)
67      nfft(3) = control_ngrid(3)
68
69      unita(1,1) = control_unita(1,1)
70      unita(2,1) = control_unita(2,1)
71      unita(3,1) = control_unita(3,1)
72      unita(1,2) = control_unita(1,2)
73      unita(2,2) = control_unita(2,2)
74      unita(3,2) = control_unita(3,2)
75      unita(1,3) = control_unita(1,3)
76      unita(2,3) = control_unita(2,3)
77      unita(3,3) = control_unita(3,3)
78
79*     **** open ELCIN binary file ****
80      if (taskid.eq.MASTER) then
81         filename = control_output_psi()
82         call util_file_name_noprefix(filename,.false.,
83     >                                .false.,
84     >                       full_filename)
85         l = index(full_filename,' ') -1
86         if (control_print(print_medium)) then
87             write(luout,1210) full_filename(1:l)
88         end if
89 1210    FORMAT(/' output psi filename:',A)
90
91         call openfile(6,full_filename,l,'w',l)
92         call iwrite(6,version,1)
93         call iwrite(6,nfft,3)
94         call dwrite(6,unita,9)
95
96         if (spin_orbit) then        !** if spin_orbit --> writeout ispin=3 **
97            ispin0 = 3
98         else
99            ispin0 = ispin
100         end if
101         call iwrite(6,ispin0,1)
102         call iwrite(6,ne,2)
103         call iwrite(6,nbrillioun,1)
104         call iwrite(6,occupation,1)
105      end if
106
107*     *** write out 3d blocks ***
108      do nb=1,nbrillioun
109         call K1dB_ktoqp(nb,nbq,pk)
110         do n=1,(ne(1)+ne(2))
111         pj = 0
112         q = n
113         if ((pj.eq.taskid_j).and.(pk.eq.taskid_k)) then
114            indx = cpsi_data_get_ptr(psi2_tag,nbq,q)
115            call Cram_c_Copy(nbq,dbl_mb(indx),dcpl_mb(tmp2(1)))
116            call Cram_c_unpack(nbq,dcpl_mb(tmp2(1)))
117         end if
118         call C3dB_c_write(1,6,dcpl_mb(tmp2(1)),
119     >                         dcpl_mb(tmp(1)),pj,pk)
120         end do
121      end do
122
123*     **** write the occupations - needs to be reworked for parallel nbrill,ne ****
124      if (occupation.gt.0) then
125         l = cpsi_data_get_next(psi2_tag)
126         do nb=1,nbrillioun
127            call K1dB_ktoqp(nb,nbq,pk)
128            call dcopy((ne(1)+ne(2)),0.0d0,0,dbl_mb(tmpe(1)),1)
129            if (pk.eq.taskid_k) then
130               call dcopy((ne(1)+ne(2)),
131     >                 dbl_mb(cpsi_data_get_chnk(l,nbq)),1,
132     >                 dbl_mb(tmpe(1)),1)
133            end if
134            call K1dB_Vector_SumAll((ne(1)+ne(2)),dbl_mb(tmpe(1)))
135            if (taskid.eq.MASTER) then
136              call dwrite(6,dbl_mb(tmpe(1)),(ne(1)+ne(2)))
137            end if
138         end do
139      end if
140
141
142*     *** close ELCIN binary file ***
143      if (taskid.eq.MASTER) then
144        call closefile(6)
145      end if
146
147c*     **** write the brillioun zone vectors ****
148c      if (taskid.eq.MASTER) then
149c         l = index(filename,' ') -1
150c         filename = filename(1:l)//"_kvectors"
151c         call util_file_name_noprefix(filename,.false.,
152c     >                                .false.,
153c     >                       full_filename)
154c         l = index(full_filename,' ') -1
155c         call openfile(6,full_filename,l,'w',l)
156c         call iwrite(6,nbrillioun,1)
157c      end if
158c      do nb=1,nbrillioun
159c         f0 = brillioun_weight_brdcst(nb)
160c         f1 = brillioun_ks_brdcst(1,nb)
161c         f2 = brillioun_ks_brdcst(2,nb)
162c         f3 = brillioun_ks_brdcst(3,nb)
163c         f4 = brillioun_k_brdcst(1,nb)
164c         f5 = brillioun_k_brdcst(2,nb)
165c         f6 = brillioun_k_brdcst(3,nb)
166c         if (taskid.eq.MASTER) then
167c         end if
168c      end do
169c
170c*     *** close ELCIN.brill binary file ***
171c      if (taskid.eq.MASTER) then
172c        call closefile(6)
173c      end if
174
175      value =           BA_pop_stack(tmpe(2))
176      value = value.and.BA_pop_stack(tmp2(2))
177      value = value.and.BA_pop_stack(tmp(2))
178      if (.not. value) call errquit('error popping stack',0, MA_ERR)
179      return
180      end
181
182
183*
184
185*     ***********************************
186*     *             			*
187*     *        cpsi_write_filename	*
188*     *             			*
189*     ***********************************
190
191      subroutine cpsi_write_filename(filename,
192     >                spin_orbit,ispin,ne,nbrillioun,psi2_tag)
193      implicit none
194      character*(*) filename
195      logical spin_orbit
196      integer ispin,ne(2),nbrillioun
197      integer psi2_tag
198
199#include "bafdecls.fh"
200#include "errquit.fh"
201
202
203*    *** local variables ***
204      integer version,l,ispin0
205      integer nfft3d
206      integer nfft(3)
207      real*8  unita(3,3)
208      character*255 full_filename
209
210      integer MASTER,taskid,taskid_j,taskid_k
211      parameter(MASTER=0)
212      integer n,nb,indx
213      integer q,pj,nbq,pk
214
215c     complex*16 tmp(*)
216      integer tmp(2),tmp2(2),occupation
217      logical value
218
219*     ***** local functions ****
220      double precision control_unita
221      external         control_unita
222      integer  control_ngrid,cpsi_data_get_ptr,cpsi_data_get_next
223      external control_ngrid,cpsi_data_get_ptr,cpsi_data_get_next
224
225      call Parallel_taskid(taskid)
226      call Parallel3d_taskid_j(taskid_j)
227      call Parallel3d_taskid_k(taskid_k)
228      call C3dB_nfft3d(1,nfft3d)
229      occupation = cpsi_data_get_next(psi2_tag)
230
231      value = BA_push_get(mt_dcpl,nfft3d,'tmp',tmp(2),tmp(1))
232      if (.not. value) call errquit('out of stack memory',0, MA_ERR)
233      value = BA_push_get(mt_dcpl,nfft3d,'tmp2',tmp2(2),tmp2(1))
234      if (.not. value) call errquit('out of stack memory',0, MA_ERR)
235
236      version = 5 ! version 5 - band, version 3 and 4 - pspw
237      nfft(1) = control_ngrid(1)
238      nfft(2) = control_ngrid(2)
239      nfft(3) = control_ngrid(3)
240
241      unita(1,1) = control_unita(1,1)
242      unita(2,1) = control_unita(2,1)
243      unita(3,1) = control_unita(3,1)
244      unita(1,2) = control_unita(1,2)
245      unita(2,2) = control_unita(2,2)
246      unita(3,2) = control_unita(3,2)
247      unita(1,3) = control_unita(1,3)
248      unita(2,3) = control_unita(2,3)
249      unita(3,3) = control_unita(3,3)
250
251*     **** open ELCIN binary file ****
252      if (taskid.eq.MASTER) then
253         call util_file_name_noprefix(filename,.false.,
254     >                                .false.,
255     >                       full_filename)
256         l = index(full_filename,' ') -1
257         write(*,1210) full_filename(1:l)
258 1210    FORMAT(/' output psi filename:',A)
259
260         call openfile(6,full_filename,l,'w',l)
261         call iwrite(6,version,1)
262         call iwrite(6,nfft,3)
263         call dwrite(6,unita,9)
264
265         if (spin_orbit) then        !** if spin_orbit --> writeout ispin=3 **
266            ispin0 = 3
267         else
268            ispin0 = ispin
269         end if
270         call iwrite(6,ispin0,1)
271         call iwrite(6,ne,2)
272         call iwrite(6,nbrillioun,1)
273         call iwrite(6,occupation,1)
274      end if
275
276*     *** write out 3d blocks ***
277      do nb=1,nbrillioun
278         call K1dB_ktoqp(nb,nbq,pk)
279         do n=1,(ne(1)+ne(2))
280            pj=0
281            q=n
282            if ((pj.eq.taskid_j).and.(pk.eq.taskid_k)) then
283               indx = cpsi_data_get_ptr(psi2_tag,nbq,q)
284               call Cram_c_Copy(nbq,dbl_mb(indx),dcpl_mb(tmp2(1)))
285               call Cram_c_unpack(nbq,dcpl_mb(tmp2(1)))
286            endif
287            call C3dB_c_write(1,6,dcpl_mb(tmp2(1)),
288     >                         dcpl_mb(tmp(1)),pj,pk)
289         end do
290      end do
291
292*     *** close ELCIN binary file ***
293      if (taskid.eq.MASTER) then
294        call closefile(6)
295      end if
296
297      value =           BA_pop_stack(tmp2(2))
298      value = value.and.BA_pop_stack(tmp(2))
299      if (.not. value) call errquit('error popping stack',0, MA_ERR)
300      return
301      end
302
303
304*     ***************************
305*     *             		*
306*     *        cpsi_write2	*
307*     *             		*
308*     ***************************
309
310      subroutine cpsi_write2(spin_orbit,ispin,ne,nbrillioun,psi2_tag)
311      implicit none
312      logical spin_orbit
313      integer ispin,ne(2),nbrillioun
314      integer psi2_tag
315
316#include "bafdecls.fh"
317#include "errquit.fh"
318
319*    *** local variables ***
320      integer version,l,ispin0,l2
321      integer nfft3d
322      integer nfft(3)
323      real*8  unita(3,3)
324      character*50 filename
325      character*255 full_filename
326      character*255 bfull_filename
327      integer MASTER,taskid,taskid_j,taskid_k
328      parameter(MASTER=0)
329      integer n,nb,indx
330      integer q,pj,nbq,pk
331
332c     complex*16 tmp(*)
333      integer tmp(2),tmp2(2),occupation
334      logical value
335
336*     ***** local functions ****
337      character*50 control_output_psi
338      external     control_output_psi
339      double precision control_unita
340      external         control_unita
341      integer  control_ngrid,cpsi_data_get_ptr,cpsi_data_get_next
342      external control_ngrid,cpsi_data_get_ptr,cpsi_data_get_next
343
344      call Parallel_taskid(taskid)
345      call Parallel3d_taskid_j(taskid_j)
346      call Parallel3d_taskid_k(taskid_k)
347      call C3dB_nfft3d(1,nfft3d)
348      occupation = cpsi_data_get_next(psi2_tag)
349
350      value = BA_push_get(mt_dcpl,nfft3d,'tmp',tmp(2),tmp(1))
351      if (.not. value) call errquit('out of stack memory',0, MA_ERR)
352      value = BA_push_get(mt_dcpl,nfft3d,'tmp2',tmp2(2),tmp2(1))
353      if (.not. value) call errquit('out of stack memory',0, MA_ERR)
354
355      version = 5 ! version 5 - band, version 3 and 4 - pspw
356      nfft(1) = control_ngrid(1)
357      nfft(2) = control_ngrid(2)
358      nfft(3) = control_ngrid(3)
359
360      unita(1,1) = control_unita(1,1)
361      unita(2,1) = control_unita(2,1)
362      unita(3,1) = control_unita(3,1)
363      unita(1,2) = control_unita(1,2)
364      unita(2,2) = control_unita(2,2)
365      unita(3,2) = control_unita(3,2)
366      unita(1,3) = control_unita(1,3)
367      unita(2,3) = control_unita(2,3)
368      unita(3,3) = control_unita(3,3)
369
370*     **** open ELCIN binary file ****
371      if (taskid.eq.MASTER) then
372         filename = control_output_psi()
373         call util_file_name_noprefix(filename,.false.,
374     >                                .false.,
375     >                       full_filename)
376         l = index(full_filename,' ') -1
377         write(*,1210) full_filename(1:l)
378         call util_file_prefix('bmovecs',bfull_filename)
379         l2 = index(bfull_filename,' ')-1
380         write(*,1210) bfull_filename(1:l2)
381 1210    FORMAT(/' output psi filename:',A)
382
383         call openfile(6,full_filename,l,'w',l)
384         call openfile(37,bfull_filename,l2,'w',l2)
385         call iwrite(6,version,1)
386         call iwrite(6,nfft,3)
387         call dwrite(6,unita,9)
388
389         if (spin_orbit) then        !** if spin_orbit --> writeout ispin=3 **
390            ispin0 = 3
391         else
392            ispin0 = ispin
393         end if
394         call iwrite(6,ispin0,1)
395         call iwrite(6,ne,2)
396         call iwrite(6,nbrillioun,1)
397         call iwrite(6,occupation,1)
398         call iwrite(37,version,1)
399         call iwrite(37,nfft,3)
400         call dwrite(37,unita,9)
401         call iwrite(37,ispin0,1)
402         call iwrite(37,ne,2)
403         call iwrite(37,nbrillioun,1)
404         call iwrite(37,occupation,1)
405      end if
406
407*     *** write out 3d blocks ***
408      do nb=1,nbrillioun
409         call K1dB_ktoqp(nb,nbq,pk)
410         do n=1,(ne(1)+ne(2))
411            pj=0
412            q=n
413            if ((pj.eq.taskid_j).and.(pk.eq.taskid_k)) then
414               indx = cpsi_data_get_ptr(psi2_tag,nbq,q)
415               call Cram_c_Copy(nbq,dbl_mb(indx),dcpl_mb(tmp2(1)))
416               call Cram_c_unpack(nbq,dcpl_mb(tmp2(1)))
417            end if
418            call C3dB_c_write(1,6,dcpl_mb(tmp2(1)),
419     >                         dcpl_mb(tmp(1)),pj,pk)
420            call C3dB_c_write(1,37,dcpl_mb(tmp2(1)),
421     >                          dcpl_mb(tmp(1)),pj,pk)
422         end do
423      end do
424
425*     *** close ELCIN binary file ***
426      if (taskid.eq.MASTER) then
427        call closefile(6)
428        call closefile(37)
429      end if
430      value =           BA_pop_stack(tmp2(2))
431      value = value.and.BA_pop_stack(tmp(2))
432      if (.not. value) call errquit('error popping stack',0, MA_ERR)
433      return
434      end
435
436
437
438
439*     ***************************
440*     *             		*
441*     *        ecpsi_write	*
442*     *             		*
443*     ***************************
444
445      subroutine ecpsi_write(spin_orbit,ispin,ne,nbrillioun,psi2_tag)
446      implicit none
447      logical spin_orbit
448      integer ispin,ne(2),nbrillioun
449      integer psi2_tag
450
451#include "bafdecls.fh"
452#include "errquit.fh"
453#include "stdio.fh"
454#include "util.fh"
455
456*    *** local variables ***
457      integer version,l,ispin0
458      integer nfft3d
459      integer nfft(3)
460      real*8  unita(3,3)
461      character*50 filename
462      character*255 full_filename
463
464      integer MASTER,taskid,taskid_j,taskid_k
465      parameter(MASTER=0)
466      integer n,nb,indx,occupation
467      integer q,pj,nbq,pk
468
469c     complex*16 tmp(*)
470      integer tmp(2),tmp2(2),tmpe(2)
471      logical value
472
473*     ***** local functions ****
474      logical      control_print
475      external     control_print
476      character*50 control_output_epsi
477      external     control_output_epsi
478      double precision control_unita
479      external         control_unita
480      integer  control_ngrid,cpsi_data_get_chnk,cpsi_data_get_next
481      integer  cpsi_data_get_ptr
482      external control_ngrid,cpsi_data_get_chnk,cpsi_data_get_next
483      external cpsi_data_get_ptr
484
485      call Parallel_taskid(taskid)
486      call Parallel3d_taskid_j(taskid_j)
487      call Parallel3d_taskid_k(taskid_k)
488      call C3dB_nfft3d(1,nfft3d)
489      occupation = -1
490
491      value = BA_push_get(mt_dcpl,nfft3d,'tmp',tmp(2),tmp(1))
492      if (.not. value) call errquit('out of stack memory',0, MA_ERR)
493      value = BA_push_get(mt_dcpl,nfft3d,'tmp2',tmp2(2),tmp2(1))
494      if (.not. value) call errquit('out of stack memory',0, MA_ERR)
495      value = BA_push_get(mt_dbl,(ne(1)+ne(2)),'tmpe',tmpe(2),tmpe(1))
496      if (.not. value) call errquit('out of stack memory',0, MA_ERR)
497
498      version = 5 ! version 5 - band, version 3 and 4 - pspw
499      nfft(1) = control_ngrid(1)
500      nfft(2) = control_ngrid(2)
501      nfft(3) = control_ngrid(3)
502
503      unita(1,1) = control_unita(1,1)
504      unita(2,1) = control_unita(2,1)
505      unita(3,1) = control_unita(3,1)
506      unita(1,2) = control_unita(1,2)
507      unita(2,2) = control_unita(2,2)
508      unita(3,2) = control_unita(3,2)
509      unita(1,3) = control_unita(1,3)
510      unita(2,3) = control_unita(2,3)
511      unita(3,3) = control_unita(3,3)
512
513*     **** open ELCIN binary file ****
514      if (taskid.eq.MASTER) then
515         filename = control_output_epsi()
516         call util_file_name_noprefix(filename,.false.,
517     >                                .false.,
518     >                       full_filename)
519         l = index(full_filename,' ') -1
520         if (control_print(print_medium))
521     >    write(luout,1210) full_filename(1:l)
522 1210    FORMAT(/' output psi filename:',A)
523
524         call openfile(6,full_filename,l,'w',l)
525         call iwrite(6,version,1)
526         call iwrite(6,nfft,3)
527         call dwrite(6,unita,9)
528
529         if (spin_orbit) then        !** if spin_orbit --> writeout ispin=3 **
530            ispin0 = 3
531         else
532            ispin0 = ispin
533         end if
534         call iwrite(6,ispin0,1)
535         call iwrite(6,ne,2)
536         call iwrite(6,nbrillioun,1)
537         call iwrite(6,occupation,1)
538      end if
539
540*     *** write out 3d blocks ***
541      do nb=1,nbrillioun
542         call K1dB_ktoqp(nb,nbq,pk)
543         do n=1,(ne(1)+ne(2))
544         pj = 0
545         q = n
546         if ((pj.eq.taskid_j).and.(pk.eq.taskid_k)) then
547            indx = cpsi_data_get_ptr(psi2_tag,nbq,q)
548            call Cram_c_Copy(nbq,dbl_mb(indx),dcpl_mb(tmp2(1)))
549            call Cram_c_unpack(nbq,dcpl_mb(tmp2(1)))
550         end if
551         call C3dB_c_write(1,6,dcpl_mb(tmp2(1)),
552     >                         dcpl_mb(tmp(1)),pj,pk)
553         end do
554      end do
555
556
557
558*     *** close ELCIN binary file ***
559      if (taskid.eq.MASTER) then
560        call closefile(6)
561      end if
562
563      value =           BA_pop_stack(tmpe(2))
564      value = value.and.BA_pop_stack(tmp2(2))
565      value = value.and.BA_pop_stack(tmp(2))
566      if (.not. value) call errquit('error popping stack',0, MA_ERR)
567      return
568      end
569
570