File: mixerifc.f90

package info (click to toggle)
elkcode 5.4.24-2
  • links: PTS, VCS
  • area: main
  • in suites: buster, sid
  • size: 12,840 kB
  • sloc: f90: 48,415; fortran: 22,457; perl: 965; makefile: 384; sh: 369; python: 105; ansic: 67
file content (68 lines) | stat: -rw-r--r-- 1,598 bytes parent folder | download
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68

! Copyright (C) 2008 J. K. Dewhurst, S. Sharma and C. Ambrosch-Draxl.
! This file is distributed under the terms of the GNU General Public License.
! See the file COPYING for license details.

subroutine mixerifc(mtype,n,v,dv,nwork,work)
use modmain
implicit none
! arguments
integer, intent(in) :: mtype,n
real(8), intent(inout) :: v(n)
real(8), intent(out) :: dv
integer, intent(inout) :: nwork
real(8), intent(inout) :: work(*)
! local variables
select case(mtype)
case(0)
! linear mixing
  if (nwork.le.0) then
    nwork=n
    return
  end if
  call mixlinear(iscl,amixpm(1),n,v,work,dv)
case(1)
! adaptive linear mixing
  if (nwork.le.0) then
    nwork=3*n
    return
  end if
  call mixadapt(iscl,amixpm(1),amixpm(2),n,v,work,work(n+1),work(2*n+1),dv)
case(3)
! Broyden mixing
  if (nwork.le.0) then
    nwork=(4+2*mixsdb)*n+mixsdb**2
    return
  end if
  call mixbroyden(iscl,n,mixsdb,broydpm(1),broydpm(2),v,work,work(2*n+1), &
   work(4*n+1),work((4+mixsdb)*n+1),work((4+2*mixsdb)*n+1),dv)
case default
  write(*,*)
  write(*,'("Error(mixerifc): mtype not defined : ",I8)') mtype
  write(*,*)
  stop
end select
return
end subroutine

subroutine getmixdata(mtype,mixdescr)
implicit none
! arguments
integer, intent(in) :: mtype
character(*), intent(out) :: mixdescr
select case(mtype)
case(0)
  mixdescr='Linear mixing'
case(1)
  mixdescr='Adaptive linear mixing'
case(3)
  mixdescr='Broyden mixing, J. Phys. A: Math. Gen. 17, L317 (1984)'
case default
  write(*,*)
  write(*,'("Error(getmixdata): mixtype not defined : ",I8)') mtype
  write(*,*)
  stop
end select
return
end subroutine