website/archive/binaries/mac/src/disembl/Tisean_3.0.1/source_f/randomize/perm/random.f

   1 c===========================================================================
   2 c
   3 c   This file is part of TISEAN
   4 c
   5 c   Copyright (c) 1998-2007 Rainer Hegger, Holger Kantz, Thomas Schreiber
   6 c
   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.
  11 c
  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.
  16 c
  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
  20 c
  21 c===========================================================================
  22 c   part of the TISEAN randomize package for constraint surrogates
  23 c   permutation scheme that swaps to randomly chosen data points
  24 c   this may also be used as a template for your own attempts
  25 c   author T. Schreiber (1999)
  26 c
  27 c-------------------------------------------------------------------
  28 c get permutation specific options
  29 c
  30       subroutine opts_permute()
  31       parameter(nx=100000)
  32       dimension nxclu(nx)
  33       character*80 filex
  34       common /permutecom/ mxclu, nxclu
  35
  36       call stcan('X',filex,' ')
  37       mxclu=0
  38       if(filex.eq." ") return
  39       open(10,file=filex,status="old",err=999)
  40  1    read(10,*,err=999,end=998) nn
  41       mxclu=mxclu+1
  42       nxclu(mxclu)=nn
  43       goto 1
  44  998  return
  45  999  write(istderr(),'(a)') "permute: cannot open "//filex
  46       stop
  47       end
  48
  49 c-------------------------------------------------------------------
  50 c print version information on permutation scheme
  51 c
  52       subroutine what_permute()
  53       call ptext("Permutation scheme: random pairs")
  54       end
  55
  56 c-------------------------------------------------------------------
  57 c print permutation specific usage message
  58 c
  59       subroutine usage_permute()
  60       call ptext("Permutation options: [-X xfile]")
  61       call popt("X", "list of indices excluded from permutation")
  62       end
  63
  64 c-------------------------------------------------------------------
  65 c initialise all that is needed for permutation scheme
  66 c
  67       subroutine permute_init()
  68       parameter(nx=100000)
  69       dimension x(nx)
  70       common nmax,cost,temp,cmin,rate,x
  71
  72       if(nmax.gt.nx) stop "permute: make nx larger."
  73       do 10 i=1,nmax
  74          call permute(n1,n2)
  75  10      call exch(n1,n2)
  76       end
  77
  78 c-------------------------------------------------------------------
  79 c find two indices n1, n2 to be exchanged, maybe using a parameter
  80 c par provided by the cooling schedule
  81 c
  82       subroutine permute(n1,n2)
  83       parameter(nx=100000)
  84       dimension nxclu(nx)
  85       common /permutecom/ mxclu, nxclu
  86       common nmax
  87       external rand
  88
  89  1    n1=min(int(rand(0.0)*nmax)+1,nmax)
  90       do 10 n=1,mxclu
  91  10      if(n1.eq.nxclu(n)) goto 1
  92  2    n2=min(int(rand(0.0)*nmax)+1,nmax)
  93       if(n2.eq.n1) goto 2
  94       do 20 n=1,mxclu
  95  20      if(n2.eq.nxclu(n)) goto 2
  96       end
  97
  98 c-------------------------------------------------------------------
  99 c given two indices n1, n2, actually perform the exchange
 100 c
 101       subroutine exch(n1,n2)
 102       parameter(nx=100000)
 103       dimension x(nx)
 104       common nmax,cost,temp,cmin,rate,x
 105
 106       h=x(n1)
 107       x(n1)=x(n2)
 108       x(n2)=h
 109       end