1 c===========================================================================
3 c This file is part of TISEAN
5 c Copyright (c) 1998-2007 Rainer Hegger, Holger Kantz, Thomas Schreiber
7 c TISEAN is free software; you can redistribute it and/or modify
8 c it under the terms of the GNU General Public License as published by
9 c the Free Software Foundation; either version 2 of the License, or
10 c (at your option) any later version.
12 c TISEAN is distributed in the hope that it will be useful,
13 c but WITHOUT ANY WARRANTY; without even the implied warranty of
14 c MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15 c GNU General Public License for more details.
17 c You should have received a copy of the GNU General Public License
18 c along with TISEAN; if not, write to the Free Software
19 c Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
21 c===========================================================================
22 c constrained randomization
23 c author T. Schreiber (1999)
24 c===========================================================================
25 parameter(nx=100000,mx=20)
27 dimension x(nx,mx), y(nx,mx), xx(nx,mx), icol(mx)
28 character*72 file, fout, comment
29 common nmax,cost,temp,cmin,rate,x
34 call whatido("constrained randomization",iverb)
38 rr=rand(ican("I",0)/real(2**22))
39 nsur=min(999,ican("n",nsur))
43 call columns(mc,mx,icol)
44 if(mcmax.eq.0) mcmax=max(1,mc)
49 isout=igetout(fout,iverb)
51 call nthstring(1,file)
52 call xreadfile(nmax,mcmax,nx,xx,nexcl,icol,file,iverb)
53 call cost_transform(nmax,mcmax,nx,xx)
54 if(file.eq."-") file="stdin"
55 if(isout.eq.1) call addsuff(fout,file,"_rnd")
56 if(nsur.gt.1) call suffix(fout,"_000")
65 if(nsur.gt.1) write(fout(index(fout," ")-3:72),'(i3.3)') isur
69 call cost_full(iv_vcost(iverb))
74 cmax=cost-temp*log(rand(0.0)) ! maximal acceptable cost
75 call cost_update(n1,n2,cmax,iaccept,iv_vmatch(iverb))
76 tnew=cool(iaccept,iend,iv_cool(iverb))
77 if(tnew.ne.temp.or.cost.lt.cmin*wr) then
79 call cost_full(iv_vcost(iverb))
80 if(iv_match(iverb).eq.1) write(istderr(),*)
81 . "cost function mismatch: ", abs((cc-cost)/cost)
84 if(cost.lt.cmin*wr) then
85 if(iv_cost(iverb).eq.1) write(istderr(),*)
86 . "after ",real(time)," steps at T=",temp," cost: ",cost
88 call cost_inverse(nmax,mcmax,nx,x,y)
89 write(comment,'(8h# cost: ,g15.5)') cost
90 call xwritecfile(nmax,mcmax,nx,y,fout,iverb,comment)
93 write(comment,'(8h# cost: ,g15.5)') cost
94 call writecfile(nmax,mcmax,nx,y,fout,iverb,comment)
101 call whatineed("[-n# -u# -I# -o outfile -l# -x# -c# -V# -h]"//
102 . " [cost opt.] [cooling opt.] [permutation opt.] file")
103 call popt("n","number of surrogates (1)")
104 call popt("u","improvement factor before write (0.9)")
105 call popt("I","seed for random numbers (0)")
106 call popt("l","number of values to be read (all)")
107 call popt("x","number of values to be skipped (0)")
108 call popt("c","column to be read (1 or file,#)")
109 call pout("file_rnd(_nnn)")
111 call ptext("Verbosity levels (add what you want):")
112 call ptext(" 1 = input/output" )
113 call ptext(" 2 = current cost if improved")
114 call ptext(" 4 = cost mismatch")
115 call ptext(" 8 = temperature etc. at cooling")
116 call ptext(" 16 = verbose cost if improved")
117 call ptext(" 32 = verbose cost mismatch")
118 write(istderr(),'()')
120 write(istderr(),'()')
122 write(istderr(),'()')
124 write(istderr(),'()')