435 lines
14 KiB
Fortran
435 lines
14 KiB
Fortran
subroutine mknet(laystr,weistr,neupop,nlay)
|
|
implicit none
|
|
! Initialize layer structure matrix exclusively.
|
|
!
|
|
! nlay: number of given layers
|
|
! neupop: neuronal population vector: contains number of neurons for each layer
|
|
! laystr: layer structure matrix
|
|
! weistr: weight structure matrix
|
|
|
|
include 'nnparams.incl'
|
|
include 'nndbg.incl'
|
|
|
|
integer laystr(3,maxlay),weistr(2,maxlay,2), neupop(maxlay)
|
|
integer nlay
|
|
|
|
integer neu_tot,wei_tot
|
|
integer n
|
|
|
|
if (dbg) write(6,'(A)') '#INITIALIZING NN....'
|
|
|
|
neu_tot=0 !number of reserved neurons spaces in memory so far
|
|
wei_tot=0 !number of reserved weight matrix element spaces in mem. so far
|
|
|
|
do n=1,nlay-1
|
|
laystr(1,n)=neupop(n) !number of neurons in layer n
|
|
laystr(2,n)=neu_tot+1 !starting pos. for layer n
|
|
laystr(3,n)=wei_tot+1 !starting pos. for W-matrix
|
|
!from layer N to N+1
|
|
neu_tot=neu_tot+neupop(n)
|
|
wei_tot=wei_tot+neupop(n)*neupop(n+1)
|
|
enddo
|
|
laystr(1,nlay)=neupop(nlay)
|
|
laystr(2,nlay)=neu_tot+1
|
|
laystr(3,nlay)=wei_tot+1 !technically oversteps boundary
|
|
!of W-vector; practically points to first bias
|
|
if (dbg) write(6,'(A)') '#PERFORMING SANITY CHECKS..'
|
|
|
|
weistr=0
|
|
|
|
! since there are no input biases, the first entry just points at
|
|
! the last few weight matrix elements
|
|
weistr(1,1,2)=wei_tot-neupop(1)+1
|
|
weistr(2,1,2)=wei_tot
|
|
|
|
do n=1,nlay-1
|
|
weistr(1,n,1)=laystr(3,n)
|
|
weistr(2,n,1)=laystr(3,n+1)-1 !end pos. of corresponding W-matrix
|
|
weistr(1,n+1,2)=weistr(2,n,2)+1 !starting pos. of bias
|
|
weistr(2,n+1,2)=weistr(2,n,2)+neupop(n+1) !end pos. of bias
|
|
enddo
|
|
|
|
|
|
if (dbg.or.vbs) write(6,*) '#PERFORMING SANITY CHECKS..'
|
|
|
|
if (neupop(1).gt.maxnin) then
|
|
write(6,'("ERROR: NETWORK INPUT EXCEEDS MAXIMUM.")')
|
|
write(6,'(I3.3,"/",I3.3)') neupop(1),maxnin
|
|
stop 1
|
|
else if (neupop(nlay).gt.maxnout) then
|
|
write(6,'("ERROR: NETWORK OUTPUT EXCEEDS MAXIMUM.")')
|
|
write(6,'(I3.3,"/",I3.3)') neupop(nlay),maxnout
|
|
stop 1
|
|
else if (nlay.gt.maxlay) then
|
|
write(6,'("ERROR: LAYER COUNT EXCEED MAXIMUM.")')
|
|
write(6,'(I2.2,"/",I2.2)') nlay,maxlay
|
|
stop 1
|
|
endif
|
|
do n=2,nlay-1
|
|
if (neupop(1).gt.maxnin) then
|
|
write(6,'("ERROR: HIDDEN LAYER SIZE EXCEEDS MAXIMUM.")')
|
|
write(6,'(I3.3,"/",I3.3)') neupop(n),maxneu
|
|
write(6,'("PROBLEM OCCURED IN LAYER #",I2.2)') n
|
|
stop 1
|
|
endif
|
|
enddo
|
|
if (dbg) then
|
|
write(6,*) '#INIT SUCCESSFUL.'
|
|
write(6,'(A10,I6)') '#NEU_TOT =',neu_tot
|
|
write(6,'(A10,I6)') '#WEI_TOT =',wei_tot
|
|
write(6,'(A10,I6)') '#MAXWEI = ',maxwei
|
|
endif
|
|
end
|
|
|
|
!--------------------------------------------------------------------------------------
|
|
|
|
subroutine mkpars(par,spread,laystr,typop,weistr,nlay,rtype)
|
|
implicit none
|
|
! Initialize Weights and Biases exclusively.
|
|
!
|
|
!
|
|
! par: one parameter vector
|
|
! spread: factor by which each weight or bias spreads around 0
|
|
! nlay: number of given layers
|
|
! rtype: type of random distribution:
|
|
! 0 -- even distribution
|
|
! 1 -- angular uniform distribution,
|
|
! vectors normalized to sqrt(wb_end)
|
|
! 2 -- even normalized distribution:
|
|
! std. deviation of weighted sums is 1.
|
|
! Ignores act.
|
|
! 3 -- same as 2, but more exact for 2+ hidden layer
|
|
! networks
|
|
!
|
|
! Initialization of the random numbers must have been done before
|
|
! call vranf(par(1:1),0,seed,iout) is used for the initialization
|
|
|
|
|
|
include 'nnparams.incl'
|
|
include 'params.incl'
|
|
include 'common.incl'
|
|
include 'nndbg.incl'
|
|
|
|
double precision par(wbcap),spread(wbcap)
|
|
integer laystr(3,maxlay),typop(maxtypes,maxlay)
|
|
integer weistr(2,maxlay,2)
|
|
integer rtype,nlay
|
|
|
|
integer wb_end
|
|
|
|
integer k
|
|
|
|
! evalute length of entire weight- and bias-vector together
|
|
wb_end=weistr(2,nlay,2)
|
|
|
|
if (rtype.eq.0) then
|
|
! generate random numbers within [0,1)
|
|
call vranf(par,wb_end,0,iout)
|
|
! shift range accordingly
|
|
do k=1,wb_end
|
|
par(k)=(par(k)-0.5D0)*2.0D0*spread(k)
|
|
enddo
|
|
else if (rtype.eq.1) then
|
|
! generate random numbers
|
|
call nnorm_grv(par,wb_end,1)
|
|
! shift range accordingly
|
|
do k=1,wb_end
|
|
par(k)=par(k)*2.0D0*spread(k)
|
|
enddo
|
|
else if (rtype.eq.2) then
|
|
! generate random numbers within [0,1)
|
|
call vranf(par,wb_end,0,iout)
|
|
! shift range to [-1,1)
|
|
do k=1,wb_end
|
|
par(k)=(par(k)-0.5D0)*2.0D0
|
|
enddo
|
|
call normalize_wbvec(par,laystr,weistr,nlay)
|
|
else if (rtype.eq.3) then
|
|
! generate random numbers within [0,1)
|
|
call vranf(par,wb_end,0,iout)
|
|
! shift range to [-1,1)
|
|
do k=1,wb_end
|
|
par(k)=(par(k)-0.5D0)*2.0D0
|
|
enddo
|
|
call normalize_wbvec_nonlin(par,laystr,weistr,typop,nlay)
|
|
par=par*spread
|
|
else
|
|
stop 'ERROR: INVALID TYPE OF RANDOM NUMBER GENERATION (RTYPE)'
|
|
endif
|
|
|
|
end
|
|
|
|
!--------------------------------------------------------------------------------------
|
|
subroutine normalize_wbvec_nonlin(par,laystr,weistr,typop,nlay)
|
|
implicit none
|
|
! Normalizes weights and biases in such a way that the standard
|
|
! deviation of a layer's weighted sums is 1, assuming par has been
|
|
! initialized with evenly distributed values from -1 to 1.
|
|
!
|
|
! For ANNs with 1 hidden layer, this routine is exact. For
|
|
! multilayered ANNs, this routine approximates distortions from ANN
|
|
! activation functions linearly. The "radius" of the even
|
|
! distribution is given for a layer k by:
|
|
!
|
|
! sqrt(3)/sqrt(1+N_eff(k-1)),
|
|
!
|
|
! where N_eff(k-1) is the effective number of neurons of the
|
|
! previous layer. It is given by the sum over all squared
|
|
! activation function derivatives at 0 of the previous layer.
|
|
! Hence, for f(x)=x for all neurons in a given layer N_eff is the
|
|
! actual number of neurons in that layer.
|
|
!
|
|
include 'nnparams.incl'
|
|
include 'params.incl'
|
|
include 'nncommon.incl'
|
|
include 'common.incl'
|
|
include 'nndbg.incl'
|
|
|
|
double precision par(wbcap)
|
|
integer laystr(3,maxlay),typop(maxtypes,maxlay)
|
|
integer weistr(2,maxlay,2)
|
|
integer nlay
|
|
|
|
double precision L(maxneu),deriv(maxneu),wbfact(maxlay)
|
|
double precision effnum
|
|
integer neu_out
|
|
|
|
integer j,n
|
|
|
|
! input layer has a trivial activation by def.
|
|
wbfact(1)=dsqrt(3.0d0)/dsqrt(dble(laystr(1,1)+1))
|
|
|
|
do n=2,nlay
|
|
L=0
|
|
deriv=0
|
|
neu_out=laystr(1,n)
|
|
|
|
! sum over squared derivatives
|
|
call neurons(neu_out,L,deriv,typop(1,n))
|
|
effnum=0.0d0
|
|
do j=1,neu_out
|
|
effnum=effnum+deriv(j)**2
|
|
enddo
|
|
|
|
wbfact(n)=dsqrt(3.0d0)/dsqrt(effnum+1.0d0)
|
|
enddo
|
|
|
|
do n=1,nlay-1
|
|
do j=weistr(1,n,1),weistr(2,n,1)
|
|
par(j)=par(j)*wbfact(n)
|
|
enddo
|
|
enddo
|
|
do n=2,nlay
|
|
do j=weistr(1,n,2),weistr(2,n,2)
|
|
par(j)=par(j)*wbfact(n-1)
|
|
enddo
|
|
enddo
|
|
|
|
end
|
|
|
|
|
|
!--------------------------------------------------------------------------------------
|
|
|
|
subroutine normalize_wbvec(par,laystr,weistr,nlay)
|
|
implicit none
|
|
! Normalizes weights and biases in such a way that the standard
|
|
! deviation of a layer's weighted sums is 1, assuming par has been
|
|
! initialized with evenly distributed values from -1 to 1.
|
|
!
|
|
! For ANNs with 1 hidden layer, this routine is exact. For
|
|
! multilayered ANNs, this routine neglects distortion effects of the
|
|
! activation functions.
|
|
include 'nnparams.incl'
|
|
include 'params.incl'
|
|
include 'nncommon.incl'
|
|
include 'common.incl'
|
|
include 'nndbg.incl'
|
|
|
|
double precision par(wbcap)
|
|
integer laystr(3,maxlay)
|
|
integer weistr(2,maxlay,2)
|
|
integer nlay
|
|
|
|
integer j,n
|
|
|
|
do n=1,nlay-1
|
|
do j=weistr(1,n,1),weistr(2,n,1)
|
|
par(j)=par(j)*dsqrt(3.0d0)/dsqrt(dble(laystr(1,n)+1))
|
|
enddo
|
|
enddo
|
|
do n=2,nlay
|
|
do j=weistr(1,n,2),weistr(2,n,2)
|
|
par(j)=par(j)*dsqrt(3.0d0)
|
|
> /dsqrt(dble(laystr(1,n-1)+1))
|
|
enddo
|
|
enddo
|
|
|
|
end
|
|
|
|
!--------------------------------------------------------------------------------------
|
|
|
|
subroutine mkmodpars(par,act,spread,nlay,thresh)
|
|
implicit none
|
|
! Assume weights and biases to be mostly initialized and only change
|
|
! those that have values that which have an absolute value less than
|
|
! thresh. Inactive weights and biases remain untouched.
|
|
!
|
|
!
|
|
! par: one parameter vector
|
|
! spread: factor by which each weight or bias spreads around 0
|
|
! act: activities for the parameter vector par. Inactive
|
|
! parameters remain untouched
|
|
! nlay: number of given layers
|
|
! thresh: minimum absolute value for a parameter to be
|
|
! overwritten. If negative, default to 0.0d0.
|
|
! Initialization of the random numbers must have been done before
|
|
! call vranf(par(1:1),0,seed,iout) is used for the initialization
|
|
|
|
|
|
include 'nnparams.incl'
|
|
include 'params.incl'
|
|
include 'common.incl'
|
|
include 'nndbg.incl'
|
|
|
|
double precision par(wbcap),spread(wbcap),thresh
|
|
integer act(wbcap)
|
|
integer nlay
|
|
|
|
double precision ranval(2)
|
|
double precision minval
|
|
integer wb_end
|
|
|
|
integer k
|
|
|
|
if (thresh.le.0.0d0) then
|
|
minval=0.0d0
|
|
else
|
|
minval=thresh
|
|
endif
|
|
|
|
! evalute length of entire weight- and bias-vector together
|
|
wb_end=pst(1,2*(nlay-1))+pst(2,2*(nlay-1))-1
|
|
|
|
do k=1,wb_end
|
|
! check whether par(k) needs changing.
|
|
if ((act(k).ne.0).and.(dabs(par(k)).le.minval)) then
|
|
! generate random number within [0,1)
|
|
call vranf(ranval,1,0,iout)
|
|
! shift range accordingly
|
|
par(k)=(ranval(1)-0.5D0)*2.0D0*spread(k)
|
|
endif
|
|
enddo
|
|
|
|
end
|
|
!--------------------------------------------------------------------------------------
|
|
|
|
subroutine mkrefset_tail(pat_in,pat_out,ref_in,ref_out,npat,nref)
|
|
implicit none
|
|
! Split off the nref last pattern pairs and make them reference data.
|
|
!
|
|
! npat: number of given pattern pairs. Will be reduced by nref.
|
|
! nref: number of reference patterns
|
|
!
|
|
! pat_in: input patterns
|
|
! pat_out: desired output patterns
|
|
! pat_*(i,N): value of ith in-/output neuron for pattern N
|
|
!
|
|
! ref_*: in analogy to pat_in/out for reference data (convergence tests)
|
|
|
|
include 'nnparams.incl'
|
|
include 'nncommon.incl'
|
|
include 'nndbg.incl'
|
|
|
|
double precision pat_in(maxnin,maxpats),pat_out(maxpout,maxpats)
|
|
double precision ref_in(maxnin,maxpats),ref_out(maxpout,maxpats)
|
|
integer npat,nref
|
|
|
|
integer j,k
|
|
|
|
if (nref.ge.npat) then
|
|
write(6,*) 'ERROR: MKREFSET: NREF TOO LARGE FOR GIVEN NPAT'
|
|
stop 1
|
|
else if (nref.lt.0) then
|
|
write(6,*) 'ERROR: MKREFSET: NEGATIVE NREF'
|
|
stop 1
|
|
endif
|
|
|
|
! update npat
|
|
npat=npat-nref
|
|
|
|
! Copy tail to ref-vectors
|
|
do j=1,nref
|
|
do k=1,maxnin
|
|
ref_in(k,j)=pat_in(k,npat+j)
|
|
enddo
|
|
do k=1,maxpout
|
|
ref_out(k,j)=pat_out(k,npat+j)
|
|
enddo
|
|
enddo
|
|
end
|
|
|
|
!--------------------------------------------------------------------------------------
|
|
|
|
subroutine mkweights(wterr,npat,pat_out)
|
|
implicit none
|
|
! Rescale weights depending on their kind (energy or state
|
|
! composition) based on a simplified form of the hybrid
|
|
! diabatisation weighting scheme. System specific rescaling is done
|
|
! by calling nnweights() for each pattern weight vector.
|
|
!
|
|
! npat: number of given pattern pairs
|
|
! wterr: weight factors for each element of the error vector e
|
|
!
|
|
! vecnum: number of state vectors in model
|
|
|
|
include 'params.incl'
|
|
include 'common.incl'
|
|
include 'nnparams.incl'
|
|
include 'nncommon.incl'
|
|
include 'nndbg.incl'
|
|
|
|
integer npat
|
|
double precision wterr(maxpout,npat),pat_out(maxpout,npat)
|
|
|
|
character*32 fname
|
|
integer j,k,n,total
|
|
|
|
fname = trim(nnfdir) // 'weights.out'
|
|
|
|
write(6,stdfmt) 'Applying weighting scheme..'
|
|
|
|
open(nnunit,file=trim(fname),status='replace')
|
|
write(nnunit,stdfmt) '# This is the established weighting scheme'
|
|
> // ' for the present fit.'
|
|
write(nnunit,stdfmt) '# Weights appear in the order given by the'
|
|
> // 'original input.'
|
|
write(nnunit,stdfmt) '# Weights given here are *not* normalized.'
|
|
|
|
write(nnunit,newline)
|
|
|
|
do k=1,npat
|
|
!call nnweight(wterr(1,k),pat_out(1,k)) ! JP no need for the weighting schemme used for energy
|
|
! right now wter is set to 1.0 in nnweights
|
|
call nnweight(wterr(1,k))
|
|
enddo
|
|
pat_out = pat_out ! to avoid unsed complain durinng compilation
|
|
total=0
|
|
do k=1,sets
|
|
write(nnunit,'(A10,I6.5)') '# Scan Nr.',k
|
|
do j=1,ndata(k)
|
|
total=total+1
|
|
do n=1,inp_out
|
|
write(nnunit,'(ES25.15)',advance='NO')
|
|
> wterr(n,total)
|
|
enddo
|
|
write(nnunit,newline)
|
|
enddo
|
|
enddo
|
|
|
|
close(nnunit)
|
|
write(6,'(A)') 'Error weight generation successful.'
|
|
write(6,'(A)') 'Wrote (not normalized) weights to '''
|
|
> // trim(fname) // '''.'
|
|
|
|
end
|