From 7f68430cb663e769a73d4e34ca22aebd2fd78812 Mon Sep 17 00:00:00 2001
From: Dave Whipp <dwhipp@dal.ca>
Date: Thu, 5 May 2011 14:30:13 +0000
Subject: [PATCH] Added sp11 and sp12 as input options for the defining the
surface geometry
---
src/create_surfaces.f90 | 13 ++++++++-----
src/module_definitions.f90 | 6 +++---
src/read_input_file.f90 | 8 ++++++++
3 files changed, 19 insertions(+), 8 deletions(-)
diff --git a/src/create_surfaces.f90 b/src/create_surfaces.f90
index 471470bb..f8b49777 100644
--- a/src/create_surfaces.f90
+++ b/src/create_surfaces.f90
@@ -50,7 +50,7 @@ integer,dimension(:), allocatable :: hulltriangles,vis_tlist,vis_elist
integer,dimension(:), allocatable :: add_tlist,nb
integer,dimension(:,:),allocatable :: vertices,neighbour,nn
logical*1,dimension(:),allocatable :: lt_work,ln_work
-double precision :: sp01,sp02,sp03,sp04,sp05,sp06,sp07,sp08,sp09,sp10
+double precision :: sp01,sp02,sp03,sp04,sp05,sp06,sp07,sp08,sp09,sp10,sp11,sp12
double precision :: x1,x2,x3,y1,y2,y3,z1,z2,z3,xne,yne,zne,xyzn
double precision :: delta,eps,epsil,xtemp,ytemp,delta_angle,delta_radius,angle
double precision :: xx,yy,xini,xend,yini,yend,dx,dy,dist
@@ -133,6 +133,8 @@ else
sp08=surface%sp08
sp09=surface%sp09
sp10=surface%sp10
+ sp11=surface%sp11
+ sp12=surface%sp12
surface_type=surface%surface_type
levelt=surface%levelt
@@ -430,7 +432,7 @@ else
!----------------------------------------------|
call zpoints (surface%nsurface,surface%x,surface%y,surface%z,surface_type, &
- sp01,sp02,sp03,sp04,sp05,sp06,sp07,sp08,sp09,sp10)
+ sp01,sp02,sp03,sp04,sp05,sp06,sp07,sp08,sp09,sp10,sp11,sp12)
!----------------------------------------------|
!----computing the normals---------------------|
@@ -452,7 +454,7 @@ end subroutine create_surf
!------------------------------------------------------------------------------|
subroutine zpoints (ns,x,y,z,surface_type, &
- sp01,sp02,sp03,sp04,sp05,sp06,sp07,sp08,sp09,sp10)
+ sp01,sp02,sp03,sp04,sp05,sp06,sp07,sp08,sp09,sp10,sp11,sp12)
!------------------------------------------------------------------------------|
!(((((((((((((((( Purpose of the routine ))))))))))))))))))))))))))))))))))))))
@@ -462,7 +464,7 @@ subroutine zpoints (ns,x,y,z,surface_type, &
! ns is the number of points
! x,y,z are the coordinates of the points
! surface_type is the type of surface under consideration
-! sp01..sp10 are the surface parameters (not all are necessarily used)
+! sp01..sp12 are the surface parameters (not all are necessarily used)
!------------------------------------------------------------------------------|
!(((((((((((((((( declaration of the subroutine arguments ))))))))))))))))))))
!------------------------------------------------------------------------------|
@@ -474,7 +476,7 @@ implicit none
integer ns
double precision x(ns),y(ns),z(ns)
integer surface_type
-double precision sp01,sp02,sp03,sp04,sp05,sp06,sp07,sp08,sp09,sp10
+double precision sp01,sp02,sp03,sp04,sp05,sp06,sp07,sp08,sp09,sp10,sp11,sp12
!------------------------------------------------------------------------------|
!(((((((((((((((( declaration of the subroutine internal variables )))))))))))))
@@ -665,6 +667,7 @@ select case(surface_type)
! sp08 is the angle of the first kink (theta1)
! sp09 is the y position of the second kink (y2)
! sp10 is the angle of the second kink (theta2)
+ ! sp11 is the dip angle of the outer 'flat' region
psi=sp03*pi/180.d0
theta1=sp08*pi/180.d0
theta2=sp10*pi/180.d0
diff --git a/src/module_definitions.f90 b/src/module_definitions.f90
index 770cd457..2754db61 100644
--- a/src/module_definitions.f90
+++ b/src/module_definitions.f90
@@ -27,8 +27,8 @@ module definitions
integer nsurface,nt
double precision,dimension(:),pointer::r,s,x,y,z,xn,yn,zn,u,v,w
integer levelt,itype,material,surface_type
- double precision :: sp01,sp02,sp03,sp04,sp05,sp06,sp07,sp08,sp09,sp10
- double precision activation_time
+ double precision :: sp01,sp02,sp03,sp04,sp05,sp06,sp07,sp08,sp09,sp10,sp11
+ double precision :: sp12,activation_time
integer,dimension(:,:),pointer::icon
logical rand
integer leveloct
@@ -330,7 +330,7 @@ module definitions
end interface
-!=====[ZISODISP================================================================
+!=====[ZISODISP]================================================================
! type ziso is used to store the isostatic deflection of the basal boundary
! due to isostasy
diff --git a/src/read_input_file.f90 b/src/read_input_file.f90
index c6fe631f..242e9257 100644
--- a/src/read_input_file.f90
+++ b/src/read_input_file.f90
@@ -422,6 +422,10 @@ do i=1,params%ns
if (iproc==0) call scanfile (params%infile,'surface_param_09_'//cm(il:3),surface(i)%sp09,ires)
surface(i)%sp10 = 0.d0
if (iproc==0) call scanfile (params%infile,'surface_param_10_'//cm(il:3),surface(i)%sp10,ires)
+ surface(i)%sp11 = 0.d0
+ if (iproc==0) call scanfile (params%infile,'surface_param_11_'//cm(il:3),surface(i)%sp11,ires)
+ surface(i)%sp12 = 0.d0
+ if (iproc==0) call scanfile (params%infile,'surface_param_12_'//cm(il:3),surface(i)%sp12,ires)
call mpi_bcast(surface(i)%sp01,1,mpi_double_precision,0,mpi_comm_world,ierr)
call mpi_bcast(surface(i)%sp02,1,mpi_double_precision,0,mpi_comm_world,ierr)
@@ -433,6 +437,8 @@ do i=1,params%ns
call mpi_bcast(surface(i)%sp08,1,mpi_double_precision,0,mpi_comm_world,ierr)
call mpi_bcast(surface(i)%sp09,1,mpi_double_precision,0,mpi_comm_world,ierr)
call mpi_bcast(surface(i)%sp10,1,mpi_double_precision,0,mpi_comm_world,ierr)
+ call mpi_bcast(surface(i)%sp11,1,mpi_double_precision,0,mpi_comm_world,ierr)
+ call mpi_bcast(surface(i)%sp12,1,mpi_double_precision,0,mpi_comm_world,ierr)
write(threadinfo%Logunit,*) 'surface ',i,'sp01',surface(i)%sp01
write(threadinfo%Logunit,*) 'surface ',i,'sp02',surface(i)%sp02
@@ -444,6 +450,8 @@ do i=1,params%ns
write(threadinfo%Logunit,*) 'surface ',i,'sp08',surface(i)%sp08
write(threadinfo%Logunit,*) 'surface ',i,'sp09',surface(i)%sp09
write(threadinfo%Logunit,*) 'surface ',i,'sp10',surface(i)%sp10
+ write(threadinfo%Logunit,*) 'surface ',i,'sp11',surface(i)%sp11
+ write(threadinfo%Logunit,*) 'surface ',i,'sp12',surface(i)%sp12
surface(i)%activation_time=-1.d0
if (iproc==0) call scanfile (params%infile,'activation_time_'//cm(il:3),surface(i)%activation_time,ires)
--
GitLab