subroutine bc1(kappa,l,jdim,jmax,ndim,n,x,t,io,task,u) c General information: c Template for writing subroutines bc1 that set the left boundary c condition for the unsteady heat conduction in a bar. Follow this c template to create different boundary conditions. c Copyright 1996 Leon van Dommelen c Version 1.0 Leon van Dommelen 12/16/96 c Usage information: c Subroutine bc1 can be called in three different ways depending on c the chosen value of parameter task: c 1. task = 0: c This is the standard call. In this case, bc1 should store the c correct boundary condition for the left-most mesh point into u(0,n). c During such a call, bc1 may assume that u(j,n) is correctly c defined for all j-values except j=0 and j=jmax. In addition, if c bc1 needs earlier time planes, it is allowed to assume that c u(j,n2) is stored cyclically, with earlier time planes stored at c n2 = n-1, n-2, ..., 0, ndim, ndim-1, ..., n+1, in order of c decreasing time. Subroutine bc1 may assume that all these earlier c time planes are fully defined as long as t(n2) is less than t(n). c This also means that if there are undefined time levels in c array u when calling subroutine bc1, the user must make sure that c the corresponding times are not less than t(n). c 2. task = -1: c This is the inquiry call. Subroutine bc1 must change the value c of task to the sum of: c 1: if the set value of u(0,n) varies nonlinearly with c changes in neighboring mesh point values u(j,n) or u(j,n2); c 2: if the boundary condition is periodic, with u(0,n) c physically equal to u(jmax-1,n); c 4: if the left-most computational point is a 'mirror c point' one mesh cell outside the actual bar; c 8: if the left-most point is half a cell outside the bar, so c that the bar end is in the middle of the first two points. c All other parameters should be ignored and no boundary condition c should be set. c 3. Any other value: c If task is positive, subroutine bc1 should merely initialize c itself and exit, otherwise it should do nothing at all. c Arguments: c Avoid typos implicit none c Input: conduction constant (restrictions may apply): double precision kappa c Input: length of the bar (restrictions may apply): double precision l c Input: declared maximum mesh point index in arrays x and u: integer jdim c Input: actual maximum mesh point index: integer jmax c Input: declared maximum time level index in arrays t and u: integer ndim c Input: time level to set in array u: integer n c Input: array of the x-positions of the mesh points: double precision x(0:jdim) c Input: array of the mesh times t: double precision t(0:ndim) c Input: I/O unit of an already open output file or zero: integer io c Subroutine bc1 should only write to this unit if io is positive. c Input/Output: task to perform: integer task c See the usage information above for more information on task. c Input/Output: temperature values: double precision u(0:jdim,0:ndim) c Provided task equals zero, bc1 must set the value of u(0,n). c Subroutine bc1 is here allowed to assume that all u(j,n) are c correctly defined for 1 <= j <= jmax-1. It may also assume that c the u(j,n2) for n2 not equal to n are defined for all j, provided c only that t(n2) < t(n). c If task is nonzero, u is undefined and should be ignored. c External variables and info for compiling or changing subroutine bc1: c The following utility routines from ../../../lib/util.f were used: c cwrite(module,text,io,show) writes line "text" to I/O unit c "io" if nonzero, and to the screen if "show" is nonzero. external cwrite c fatal(module,text1,text2,text3) kills the program after a c fatal error, printing the lines "text1", "text2" and "text3". external fatal c Local variables: c *************************************************************** c ******* Declare any local variables that bc1 needs here ******* c *************************************************************** c An integer to keep track of whether bc1 has been initialized: integer init c Integer codes for the various types of boundary conditions: integer nonlin,period,mirror,staggr parameter (nonlin=1,period=2,mirror=4,staggr=8) c When creating a new implementation of the boundary condition c requiring special locations x(j) of the mesh points, give the new c boundary condition its own code (powers of 2 only). Then make sure c that ../init/initx.f creates the desired mesh point x-values when c given this code. c For a complete list of all currently used codes, please see the c source code of initx. c Constants defined to make changing precision easier: double precision zero parameter (zero=0.d0) c Data statements: data init/0/ c Executable statements: c For task = -1, return the properties of the boundary condition: if(task.eq.-1)then c ********************************************************* c *** Sum the properties of the boundary condition here *** c ********************************************************* task=? return endif c Ignore other negative tasks: if(task.lt.0)return c Initialization during the first time that bc1 is called: if(init.eq.0)then c Show which subroutine is being initialized: call cwrite('bc1', & 'Initialization of left boundary condition bc1:',io,0) c ************************************************************ c *********** Briefly describe subroutine bc1 here *********** c ************************************************************ call cwrite('bc1', & '- template/bc_1.f Copyright 1996 Leon van Dommelen.',io,0) call cwrite('bc1', & '*** ERROR: You are using the template subroutine!',io,0) c ********************************************************* c ******* Read in any parameters needed over here ********* c ********************************************************* c All done with the initialization: init=1 endif c Further ignore nonzero tasks: if(task.ne.0)return c Check the sanity of some of the arguments: if(jdim.le.0 .or. ndim.lt.0 .or. & jmax.le.0 .or. n.lt.0 .or. & jmax.gt.jdim .or. n.gt.ndim)call fatal('bc1', & 'Program error: invalid arguments issued to subroutine bc1.', & 'Fix the program code.',' ') c Set the boundary point: c **************************************************************** c ************ Compute the boundary point u(0,n) here ************ c **************************************************************** u(0,n)=zero? goto 900 c Exit: c Jump here when done: 900 return end