2D Topography Model

1

I’m working on a Fortran program to create a 2D model for near-surface topography, specifically simulating a sinusoidal topography on the second layer’s top surface. The model outputs a grid to a text file, with the first layer represented by 1 and the second layer by 2. The spacing between grid points is (nx - 1). Here is the code I’ve written so far:

program Sin_Topo_Layer

implicit none

integer :: nx   ! Number of grid points in the x direction
integer :: nz   ! Number of grid points in the y direction
real :: x, z
integer :: i, j
!real :: z(ny, nx)

! integer :: topo(ny, nx)
!real :: x_min, x_max, z_min, z_max
!real :: dx, dz

real :: amplitude    ! Amplitude of the sinusoidal wave
real :: wavelength   ! Wavelength of the sinusoidal wave

character*1,allocatable,dimension(:,:slight_smile: :: icode

! Values of parameters

nx = 100
nz = 100
amplitude = 5.0  
wavelength = 20.0  
x_min = 0.0
x_max = 10.0
z_min = 0.0
z_max = 10.0



! Compute the grid spacing
dx = (x_max - x_min) / (nx - 1)
dz = (z_max - z_min) / (nz - 1)

! Allocate the topography array

    	allocate(icode(1:nx,1:nz))
        
         ! Create topography
do j = 1, nz
    do i = 1, nx
        x = (i - 1) * wavelength / (nx - 1)
        z = (j - 1) * wavelength / (nz - 1)
        z(j, i) = amplitude * sin(2.0 * 3.141592653589793 * x / wavelength)

! Assign layers based on the sinusoidal topography
if (y > z(j, i)) then
topo(j, i) = 1 ! First layer
else
topo(j, i) = 2 ! Second layer
end if
end do
end do

! Write the topography to the output file
open(2012,file='Sin_Topo_Layer.txt',status='unknown')
	do j=1,nz
		write(2012,1000)(icode(i,j),i=1,nx)
	enddo
close(2012)
1000 format(<nx>a1)

! Generate the topography using sin fuction 
!do j = 1, nz
   ! z = z_min + (j - 1) * dz
  !  do i = 1, nx
     !   x = x_min + (i - 1) * dx
     !   topo(i, j) = sin(x) * cos(z)
    !end do
!end do

! Output the data to a text file
! do j = 1, nz
!    z = z_min + (j - 1) * dz
 !   do i = 1, nx
   !     x = x_min + (i - 1) * dx
   !     if (j <= nz / 2) then
    !        topo(i, j) = 1  ! First layer
    !    else
    !        topo(i, j) = 2  ! Second layer
   !     end if
   ! end do

! end do
! deallocate(topo)’

!open(2012,file='Sin_Topo_Layer.txt',status='unknown')
!	do j=1,nz
!		write(2012,1000)(icode(i,j),i=1,nx)
!	enddo
!close(2012)
!1000 format(<nx>a1)

end program Sin_Topo_Layer

Snipaste_2024-07-04_16-10-19

1 Like
2

! Sin_Topo_Layer.f90
!
! FUNCTIONS:
! Sin_Topo_Layer - Entry point of console application.
! Create a 2D model for near-surface topography (Sin Top Layer)
! simulate a sinusoidal topography on the second layer’s top surface.
!Below is a basic Fortran code snippet that generates a 2D grid with sinusoidal topography for the second layer

!****************************************************************************
!
! PROGRAM: Sin_Topo_Layer
!
! PURPOSE: Entry point for the console application.
!
!****************************************************************************

program Sin_Topo_Layer

implicit none

    ! Parameters
integer:: nx  ! Number of points along x-axis
integer:: nz  ! Number of points along y-axis
real, parameter :: x_min = 0.0  ! Minimum x coordinate
real, parameter :: x_max = 10.0 ! Maximum x coordinate
real, parameter :: z_min = 0.0  ! Minimum y coordinate
real, parameter :: z_max = 10.0  ! Maximum y coordinate
real :: amplitude        ! Amplitude of the sinusoidal wave
real :: wavelength       ! Wavelength of the sinusoidal wave
integer:: idepth
integer :: i, j
real :: x, z, z_top
real :: dx, dz
character*1,allocatable,dimension(:,:) :: icode

! Values of parameters
nx = 1001
nz = 501
idepth = 10  ! Depth where the second layer starts
amplitude  = 10.0  
wavelength = 6.0  

! Compute the grid spacing
dx = (x_max - x_min) / (nx - 1)
dz = (z_max - z_min) / (nz - 1)

! Allocate the array for the topography

	allocate(icode(1:nx,1:nz))

! Create topography
do j = 1, nz
    z = z_min + (j - 1) * dz
    do i = 1, nx
        x = x_min + (i - 1) * dx
        z_top = amplitude * sin(2.0 * 3.141592653589793 * x / wavelength)
        if (j > nz / 2 + z_top * (nz / amplitude / 2)) then
            icode(i, j) = '1'  ! First layer
        else if (j >= idepth) then
            icode(i, j) = '2'  ! Second layer
       else
          icode(i, j) = '3'  ! Undefined or empty space
        end if
    end do
end do


! Write the topography to the output file
open(2012,file='Sin_Topo_Layer.txt',status='unknown')
	do j=1,nz
		write(2012,1000)(icode(i,j),i=1,nx)
	enddo
close(2012)

1000 format(a1)

end program Sin_Topo_Layer

Help
little help

3

I want the result to be like this

Sin_topo_layer.TXT (259.4 KB)

4

Could you explain your question in a bit more detail? Do you need help with writing the header or with filling the matrix? What are the parameters that define the sinus shape? Note: I have not studied your code in any detail, as I did not see any concrete question ;).

5

I’m working on a Fortran program to create a 2D model for near-surface topography, specifically simulating a sinusoidal topography on the second layer’s top surface (Sin—Topo—Layer). The program outputs a grid to a text file, with the first layer represented by 1 and the second layer represented by 2. The spacing between grid points is determined by (nx—1). The current Fortran code I wrote has some mistakes, and the resulting topography shape is not accurate. To clarify the expected outcome, please refer to the attached (Sin_topo_layer.TXT), which shows the correct format.

Snipaste_2024-07-04_16-10-19
Sin_topo_layer.TXT (259.4 KB)

6

Snipaste_2024-07-08_13-45-08
Snipaste_2024-07-08_13-45-08846×630 77.3 KB

7

I ran your program and all the 1’s, 2’s and 3’s occurred on a line of their own. Not what you want if the example file is what you require. I changed the format to “(1001a1)” and I got a sine wave which is much deeper than the one in the example. But it is also nicely sinusy - not as jagged as the picture you show. This helped me understand what is going wrong: the sine’s amplitude is nz/2, not 10.

1 Like