Project

General

Profile

MagOptim.f90

Zaharko, Oksana, 05/02/2012 12:40 PM

Download (11.7 KB)

 
1
Program Optimizing_MagStructure
2
   !---- Optimizing Magnetic Structure vs Integrated Intensity for a MultiDomain Case
3
   !---- Use Modules ----!
4
   use CFML_crystallographic_symmetry, only: space_group_type, Write_SpaceGroup
5
   use CFML_string_utilities,          only: u_case, Get_LogUnit
6
   use CFML_Atom_TypeDef,              only: Atom_List_Type, Write_Atom_List, mAtom_list_Type
7
   use CFML_crystal_Metrics,           only: Crystal_Cell_Type, Write_Crystal_Cell
8
   use CFML_IO_Formats,                only: Readn_set_Xtal_Structure,err_form_mess,err_form,file_list_type
9
   use CFML_Keywords_Code_Parser,      only: NP_Max,NP_Refi,V_BCon,V_Bounds,V_Name,V_Vec,&
10
                                             Allocate_VParam,Init_RefCodes, Read_RefCodes_File, &
11
                                             Write_Info_RefCodes, Err_RefCodes, err_refcodes_mess, &
12
                                             VState_to_AtomsPar
13
   use CFML_Magnetic_Symmetry
14
   use CFML_Simulated_Annealing,       only: SimAnn_Conditions_type, state_Vector_Type, Multistate_Vector_Type, &
15
                                             err_san_mess,err_SAN, Simanneal_Gen,Set_SimAnn_Cond, &
16
                                             Set_SimAnn_StateV,Write_SimAnn_Cond, Write_SimAnn_StateV, &
17
                                             Write_SimAnn_MStateV, SimAnneal_MultiConf,Set_SimAnn_MStateV, &
18
                                             Sann_Opt_MultiConf
19
   use cost_magfunctions                 
20
   use prep_input
21
   
22
   !---- Local Variables ----!
23
   implicit none
24

    
25
   !---- List of public subroutines ----!
26

    
27
   type (SimAnn_Conditions_type),save :: c
28
   type (state_Vector_Type),save      :: vs
29
   type (Multistate_Vector_Type)      :: mvs
30
   type (Atom_list_Type)              :: A
31
   
32
   character(len=256)                 :: filcod     !Name of the input file
33
   character(len=256)                 :: filhkl     !Name of the hkl-file
34
   character(len=256)                 :: filrest    !Name of restraints file
35
   character(len=256)                 :: line       !Text line
36
   character(len=256)                 :: fst_cmd    !Commands for FP_Studio
37
   integer                            :: Num, ier,i,j,k,n, i_cfl
38
   real                               :: start,fin
39
   integer                            :: narg,iargc, n_ini,n_end
40
   Logical                            :: esta, arggiven=.false.,sthlgiven=.false., &
41
                                         fst_out=.false., local_opt=.false., rest_file=.false.
42

    
43
    !---- Arguments on the command line ----!
44
    narg=COMMAND_ARGUMENT_COUNT()
45

    
46
    if(narg > 0) then
47
            call GET_COMMAND_ARGUMENT(1,filcod)
48
            arggiven=.true.
49
    end if
50

    
51
    write (unit=*,fmt="(/,/,6(a,/))")                                                  &
52
            "            ------ PROGRAM FOR OPTIMIZING Magnetic STRUCTURE ------"            , &
53
            "                    ---- Version 1.1 February-2012 ----"                        , &
54
            "    *******************************************************************************"  , &
55
            "    * Presently optimizes a magnetic structure given in *.CFL file           "  , &
56
            "    * against magnetic intensities listed in *.INT file                      "  , &
57
            "    *******************************************************************************"  , &
58
            "                             (OZ-February-2012 )"
59
   write (unit=*,fmt=*) " "
60

    
61
   if(.not. arggiven) then
62
     write(unit=*,fmt="(a)",advance='no') " => Code of the file xx.cfl (give xx): "
63
     read(unit=*,fmt="(a)") filcod
64
     if(len_trim(filcod) == 0) stop
65
   end if
66

    
67
   open(unit=lun,file=trim(filcod)//".out", status="replace",action="write")
68
   write(unit=lun,fmt="(/,/,6(a,/))")                                                  &
69
            "            ------ PROGRAM FOR OPTIMIZING Magnetic STRUCTURE ------"            , &
70
            "                    ---- Version 1.1 February-2012 ----"                        , &
71
            "    *******************************************************************************"  , &
72
            "    * Presently optimizes a magnetic structure given in *.CFL file           "  , &
73
            "    * against magnetic intensities listed in *.INT file                      "  , &
74
            "    *******************************************************************************"  , &
75
            "                             (OZ-February-2012 )"
76

    
77
   inquire(file=trim(filcod)//".cfl",exist=esta)
78
   if( .not. esta) then
79
     write(unit=*,fmt="(a)") " File: "//trim(filcod)//".cfl doesn't exist!"
80
     stop
81
   end if
82

    
83
  
84
!!----
85
!!---- reading from CFL file
86
!!----
87
   call Readn_Set_Xtal_Structure(trim(filcod)//".cfl",Cell,SpG,A,Mode="CFL",file_list=file_cfl)
88

    
89
   If(err_form) then
90
     write(unit=*,fmt="(a)") trim(err_form_mess)
91
   Else
92
     call Write_Crystal_Cell(Cell,lun)
93
     call Write_SpaceGroup(SpG,lun)
94
     call Write_Atom_List(A,lun=lun)
95

    
96
     n_ini=1
97
     n_end=file_cfl%nlines
98
     call Readn_Set_Magnetic_Structure(file_cfl,n_ini,n_end,MGp,mA,Mag_dom=Mag_dom,Cell=Cell) !Cell is needed for Spherical components
99
       if(err_MagSym) then
100
         write(unit=*,fmt="(a)") "   "//err_MagSym_Mess
101
         stop
102
       end if
103
     call Write_Magnetic_Structure(lun,mGp,mA,Mag_dom)
104
     NP_Max= mA%natoms*17+Mag_dom%nd*2  ! Rxyz,Ixyz (also Rm, Rphi, Rth), mPhase, C1-12 + domains
105
     call Allocate_Vparam(NP_Max)
106

    
107
     !Determine the refinement codes from vary/fix instructions
108
     call Init_RefCodes(mA)
109
     call Init_RefCodes(mag_Dom)
110
     n_ini=1
111
     n_end=file_cfl%nlines
112
     call Read_RefCodes_File(file_cfl,n_ini,n_end,mA,Spg,Mag_dom)
113

    
114
     if(Err_RefCodes) then
115
       write(unit=*,fmt="(a)") trim(err_refcodes_mess)
116
     end if
117

    
118
     call Write_Info_RefCodes(FmAtom=mA,Mag_dom=Mag_dom,MGp=mGp,Iunit=lun)
119

    
120
     !--- Read data file name(s)
121
     call Readn_Set_sqMiV_Data(file_cfl)
122

    
123
     !--- Read what cost functions will be minimized
124
     call Readn_Set_CostFunctPars(file_cfl)
125
     if(Err_cost) then
126
       write(unit=*,fmt="(a)") trim(Err_Mess_cost)
127
       stop
128
     end if
129
     call Write_CostFunctPars(lun)
130
 
131
      write(unit=*,fmt="(a)") " Normal End of Routine: Read_Cfl"
132
      write(unit=*,fmt="(a)") " Results in File: "//trim(filcod)//".out"
133

    
134

    
135
     !--- Read job and subjob (Opti, f2mag) 
136
     call Readn_Set_Job(file_cfl)
137

    
138
     !--- Set up the simulated annealing conditions
139
     call Set_SimAnn_Cond(file_cfl,c)
140
      if(err_SAN) then
141
         write(unit=*,fmt="(a)") " => Error setting Simulated Annealing conditions"
142
         write(unit=*,fmt="(a)") trim(err_san_mess)
143
        stop
144
      end if
145
     call Write_SimAnn_Cond(lun,c)
146
     call cpu_time(start)
147
     if(c%num_conf == 1) then
148
        !--- Set up the Simulated Annealing State Vector
149
         call Set_SimAnn_StateV(NP_Refi,V_BCon,V_Bounds,V_Name,V_Vec,vs)
150
          if(err_SAN) then
151
             write(unit=*,fmt="(a)") " => Error setting Simulated Annealing State Vector"
152
             write(unit=*,fmt="(a)") trim(err_san_mess)
153
            stop
154
          end if
155
         call Write_SimAnn_StateV(lun,vs,"INITIAL STATE")
156
        !--- Here Simulated Annealing is performed
157
         call Simanneal_Gen(General_Cost_function,c,vs,lun)
158
         call Write_SimAnn_StateV(lun,vs,"FINAL STATE")
159

    
160
         !---Write a CFL file
161
         call Copy_mAtom_List(mA, mA_clone)
162
         call Copy_MagDom_List(Mag_dom, Mag_dom_clone)
163
         call Get_LogUnit(i_cfl)
164
         write(unit=line,fmt="(a,f12.2)") "  Best Configuration found by Simanneal_Gen,  cost=",vs%cost
165
         open(unit=i_cfl,file=trim(filcod)//"_sol.cfl",status="replace",action="write")
166
         call Write_SOL_mCFL(i_cfl,file_cfl,mA_Clone,Mag_dom_clone,line) 
167
         close(unit=i_cfl)
168

    
169
     else   ! Multi-configuration Simulated Annealing
170

    
171
         call Set_SimAnn_MStateV(NP_Refi,c%num_conf,V_BCon,V_Bounds,V_Name,V_Vec,mvs)
172
          if(err_SAN) then
173
             write(unit=*,fmt="(a)") " => Error setting Simulated Annealing State Vector"
174
             write(unit=*,fmt="(a)") trim(err_san_mess)
175
            stop
176
          end if
177
         call Write_SimAnn_MStateV(lun,mvs,"INITIAL STATE")
178
        !--- Here Simulated Annealing is performed
179
         call SAnn_Opt_MultiConf(General_Cost_function,c,mvs,lun)
180
         call Write_SimAnn_MStateV(lun,mvs,"FINAL STATE",1)
181
         !Write CFL files
182
         call Copy_mAtom_List(mA, mA_clone)
183
         call Copy_MagDom_List(Mag_dom, Mag_dom_clone)
184
         n=mvs%npar
185
         do i=1,c%num_conf
186
           line=" "
187
           write(unit=line,fmt="(a,i2.2,a)") trim(filcod)//"_sol",i,".cfl"
188
           V_vec(1:n)=mvs%state(1:n,i)
189
           call VState_to_AtomsPar(mA_clone,mode="V",MGp=Mgp,Mag_dom=Mag_dom_clone) ! V-Passing Value
190
           call Get_LogUnit(i_cfl)
191
           open(unit=i_cfl,file=trim(line),status="replace",action="write")
192
           write(unit=line,fmt="(a,i2,a,f12.2)") "  Configuration #",i," found by SimAnneal_MultiConf,  cost=",mvs%cost(i)
193
           call Write_SOL_mCFL(i_cfl,file_cfl,mA_Clone,Mag_dom_clone,line)
194
           close(unit=i_cfl)
195
         end do
196
     end if
197

    
198
     call Write_FinalCost(lun)
199
     call Write_SOL_mCFL(lun,file_cfl,mA_Clone,Mag_dom_clone,line)
200

    
201
     !--- Here results are written
202
     call Write_ObsCalc(file_cfl)
203
     
204
     write(unit=*,fmt="(a)") " Normal End of: PROGRAM FOR OPTIMIZING magnetic STRUCTURES "
205
     write(unit=*,fmt="(a)") " Results in File: "//trim(filcod)//".out"
206
     call cpu_time(fin)
207
     write(unit=*,fmt="(a,f10.2,a)")  "  CPU-Time: ", fin-start," seconds"
208
     write(unit=*,fmt="(a,f10.2,a)")  "  CPU-Time: ", (fin-start)/60.0," minutes"
209
     write(unit=lun,fmt="(/,a,f10.2,a)")  "  CPU-Time: ", fin-start," seconds"
210
   end if
211

    
212
   close(unit=lun)
213
   stop
214
End Program Optimizing_MagStructure
215

    
216
Subroutine Write_FST(fst_file,v,cost) ! Is not used here, needed for CFML_optimization_san from libcrysfml
217
   !---- Arguments ----!
218
   use CFML_String_Utilities,        only: get_logunit
219
   use CFML_Keywords_Code_Parser,    only: VState_to_AtomsPar
220
   use CFML_Atom_TypeDef,            only: atom_type, atom_list_type
221
   use CFML_crystal_Metrics,         only: Crystal_Cell_Type
222
   use CFML_crystallographic_symmetry, only: space_group_type
223

    
224
   character(len=*),     intent(in):: fst_file
225
   real,dimension(:),    intent(in):: v
226
   real,                 intent(in):: cost
227

    
228
   !----- Local variables -----!
229
   type(Atom_Type)                    :: atom
230
   type (Atom_list_Type)              :: A
231
   type (Crystal_Cell_Type)           :: Cell
232
   type (space_group_type)            :: SpG
233

    
234
   integer :: lun,i,nc, ier
235
   character(len=132)                 :: file_fst,fst_cmd
236
   character(len=30), dimension(10)   :: cmds
237

    
238
   i=index(fst_file,".fst")
239
   file_fst=fst_file(1:i+3)
240
   fst_cmd=adjustl(fst_file(i+4:))
241
   nc=0
242

    
243
   do
244
     i=index(fst_cmd,";")
245
     if(i /= 0) then
246
       nc=nc+1
247
       cmds(nc)=fst_cmd(1:i-1)
248
       fst_cmd=fst_cmd(i+1:)
249
     else
250
       nc=nc+1
251
       cmds(nc)=fst_cmd
252
       exit
253
     end if
254
   end do
255

    
256
   !Update the atom parameters
257
   call VState_to_AtomsPar(A,mode="V")
258

    
259
   call get_logunit(lun)
260
   open(unit=lun, file=trim(file_fst), status="replace",action="write",position="rewind",iostat=ier)
261
   if (ier == 0) then
262
      write(unit=lun,fmt="(a)") "TITLE  FST-file generated with Write_FST"
263
      write(unit=lun,fmt="(a)") "SPACEG "//trim(Spg%SPG_Symb)
264
      write(unit=lun,fmt="(a,3f12.5,3f8.3)") "CELL ",cell%cell, cell%ang
265
      write(unit=lun,fmt="(a)") "BOX   -0.20  1.20   -0.20  1.20    -0.20  1.20 "
266
      do i=1,A%natoms
267
         write(unit=lun,fmt="(a,a,3f12.5,tr4,a)")"Atom "//A%atom(i)%lab, A%atom(i)%chemsymb, A%atom(i)%x, A%atom(i)%AtmInfo
268
      end do
269
      do i=1,nc
270
         write(unit=lun,fmt="(a)") trim(cmds(i))
271
      end do
272
      write(unit=lun,fmt="(a)") "!"
273
      call flush(lun)
274
      !flush(lun)
275
      close(unit=lun)
276
   end if
277

    
278
   return
279
End Subroutine Write_FST