debug_pml.f90

Sfdtd Simple finite-difference time-domain

! debug_pml.f90
!
!    density plot with PML layers    
!
!    Copyright (C) 2007  Paul Panserrieu, < [e-mail] >
!
!    This program is free software: you can redistribute it and/or modify
!    it under the terms of the GNU General Public License as published by
!    the Free Software Foundation, either version 3 of the License.
! 
! last modified: 19-10-2007 12:22:10 PM CEST

MODULE debug_pml

USE pml, ONLY: pml_gitter, pml_boundary
USE plot, ONLY: toch
USE fdtd_gitter, ONLY: gitter

IMPLICIT NONE

CONTAINS

SUBROUTINE plot_pml_debug(g, b, it, n, plane, ebene) 

  TYPE(gitter), INTENT(IN)           :: g
  TYPE(pml_boundary), INTENT(IN)     :: b
  INTEGER, INTENT(IN)                :: it   
  INTEGER, INTENT(IN)                :: n   
  CHARACTER, INTENT(IN)              :: plane
  INTEGER, INTENT(IN)               :: ebene

  CALL density_pml(g, b, it, ebene, plane, 'E', n, 1, 'xxx')
  CALL density_pml(g, b, it, ebene, plane, 'E', n, 2, 'xxx')
  CALL density_pml(g, b, it, ebene, plane, 'E', n, 3, 'xxx')

  CALL density_pml(g, b, it, ebene, plane, 'H', n, 1, 'xxx')
  CALL density_pml(g, b, it, ebene, plane, 'H', n, 2, 'xxx')
  CALL density_pml(g, b, it, ebene, plane, 'H', n, 3, 'xxx')

END SUBROUTINE plot_pml_debug

SUBROUTINE density_pml(g, b, t, ne, co, cf, n, component, memo)
    
  TYPE(gitter), INTENT(IN)           :: g
  TYPE(pml_boundary), INTENT(IN)     :: b
  INTEGER, INTENT(IN)                :: t   
  INTEGER, INTENT(IN)                :: ne   
  INTEGER, INTENT(IN)                :: n  
  CHARACTER, INTENT(IN)              :: co   
  CHARACTER, INTENT(IN)              :: cf   
  INTEGER, INTENT(IN)                :: component 
  CHARACTER(len=3), INTENT(IN)       :: memo
  INTEGER                            :: i1, i2
  CHARACTER, DIMENSION(1:5)          :: x
  CHARACTER(len=5)                   :: str
  CHARACTER                          :: komp

  CALL toch(x, t)
  str = x(1)//x(2)//x(3)//x(4)//x(5)
  IF (component .EQ. 1) THEN
    komp = 'x'
  ELSEIF (component .EQ. 2) THEN
    komp = 'y'
  ELSE
    komp = 'z'
  ENDIF
  OPEN(10,FILE='pml_'//g%id//'_'//co//'_'//cf//komp//'_'//memo//str//'.plt',ACTION='WRITE')

  WRITE(10,*) '$ DATA=CONTOUR'
  WRITE(10,*) '% toplabel   = " "'
  WRITE(10,*) '% subtitle   = "Zeitschritte: ', t, '"'
  WRITE(10,*) '% interp     = 2'
  WRITE(10,*) '% contfill   = on'
  WRITE(10,*) '% meshplot   = off'
  WRITE(10,*) '% nsteps     = 50'
  WRITE(10,*) '% nx   =',  n + 1 + 2 * b%layers
  WRITE(10,*) '% ny   =',  n + 1 + 2 * b%layers
  WRITE(10,*) '% xmin = ', g%nxl-b%layers, 'xmax =', g%nxgh+b%layers
  WRITE(10,*) '% ymin = ', g%nxl-b%layers, 'ymax =', g%nxgh+b%layers

  IF ( co .EQ. 'x' ) THEN
    WRITE(10,*) '% xlabel="y"'
    WRITE(10,*) '% ylabel="z"'
    IF (ne < g%nxl) THEN
      DO i2 = g%nxl-b%layers, g%nxgh+b%layers, 1
        DO i1 = g%nxl-b%layers, g%nxgh+b%layers, 1
          IF( cf .EQ. 'E'  .OR. cf .EQ. 'e') THEN
            WRITE(10,*)   b%bas_x%E(ne, i1, i2, 2 * component - 1) &
                        + b%bas_x%E(ne, i1, i2, 2 * component)
          ELSE
            WRITE(10,*)   b%bas_x%H(ne, i1, i2, 2 * component - 1) &
                        + b%bas_x%H(ne, i1, i2, 2 * component)
          ENDIF
        ENDDO
      ENDDO
    ELSEIF (ne > g%nxyh) THEN
      DO i2 = g%nxl-b%layers, g%nxgh+b%layers, 1
        DO i1 = g%nxl-b%layers, g%nxgh+b%layers, 1
          IF( cf .EQ. 'E'  .OR. cf .EQ. 'e') THEN
            WRITE(10,*)   b%top_x%E(ne, i1, i2, 2 * component - 1) &
                        + b%top_x%E(ne, i1, i2, 2 * component)
          ELSE
            WRITE(10,*)   b%top_x%H(ne, i1, i2, 2 * component - 1) &
                        + b%top_x%H(ne, i1, i2, 2 * component)
          ENDIF
        ENDDO
      ENDDO
    ELSE
      DO i2 = g%nxl-b%layers, g%nxgh+b%layers, 1
        DO i1 = g%nxl-b%layers, g%nxgh+b%layers, 1
          IF (i1 < g%nyl) THEN      
            IF( cf .EQ. 'E'  .OR. cf .EQ. 'e') THEN
              WRITE(10,*)   b%bas_y%E(ne, i1, i2, 2 * component - 1) &
                          + b%bas_y%E(ne, i1, i2, 2 * component)
            ELSE
              WRITE(10,*)   b%bas_y%H(ne, i1, i2, 2 * component - 1) &
                          + b%bas_y%H(ne, i1, i2, 2 * component)
            ENDIF
          ELSEIF (i1 > g%nyyh) THEN
            IF( cf .EQ. 'E'  .OR. cf .EQ. 'e') THEN
              WRITE(10,*)   b%top_y%E(ne, i1, i2, 2 * component - 1) &
                          + b%top_y%E(ne, i1, i2, 2 * component)
            ELSE
              WRITE(10,*)   b%top_y%H(ne, i1, i2, 2 * component - 1) &
                          + b%top_y%H(ne, i1, i2, 2 * component)
            ENDIF
          ELSEIF (i2 < g%nzl) THEN
            IF( cf .EQ. 'E'  .OR. cf .EQ. 'e') THEN
              WRITE(10,*)   b%bas_z%E(ne, i1, i2, 2 * component - 1) &
                          + b%bas_z%E(ne, i1, i2, 2 * component)
            ELSE
              WRITE(10,*)   b%bas_z%H(ne, i1, i2, 2 * component - 1) &
                          + b%bas_z%H(ne, i1, i2, 2 * component)
            ENDIF
          ELSEIF (i2 > g%nzyh) THEN
            IF( cf .EQ. 'E'  .OR. cf .EQ. 'e') THEN
              WRITE(10,*)   b%top_z%E(ne, i1, i2, 2 * component - 1) &
                          + b%top_z%E(ne, i1, i2, 2 * component)
            ELSE
              WRITE(10,*)   b%top_z%H(ne, i1, i2, 2 * component - 1) &
                          + b%top_z%H(ne, i1, i2, 2 * component)
            ENDIF
          ELSE
            IF( cf .EQ. 'E'  .OR. cf .EQ. 'e') THEN
              WRITE(10,*) g%E(ne, i1, i2, component)
            ELSE
              WRITE(10,*) g%H(ne, i1, i2, component)
            ENDIF
          ENDIF
        ENDDO
      ENDDO
    ENDIF
  ENDIF

  IF ( co .EQ. 'y' ) THEN
    WRITE(10,*) '% xlabel="x"'
    WRITE(10,*) '% ylabel="z"'
    IF (ne < g%nyl) THEN
      DO i2 = g%nxl-b%layers, g%nxgh+b%layers, 1
        DO i1 = g%nxl-b%layers, g%nxgh+b%layers, 1
          IF (i1 < g%nxl) THEN
            IF( cf .EQ. 'E'  .OR. cf .EQ. 'e') THEN
              WRITE(10,*)   b%bas_x%E(i1, ne, i2, 2 * component - 1) &
                          + b%bas_x%E(i1, ne, i2, 2 * component)
            ELSE
              WRITE(10,*)   b%bas_x%H(i1, ne, i2, 2 * component - 1) &
                          + b%bas_x%H(i1, ne, i2, 2 * component)
            ENDIF
          ELSEIF (i1 > g%nxyh) THEN
            IF( cf .EQ. 'E'  .OR. cf .EQ. 'e') THEN
              WRITE(10,*)   b%top_x%E(i1, ne, i2, 2 * component - 1) &
                          + b%top_x%E(i1, ne, i2, 2 * component)
            ELSE
              WRITE(10,*)   b%top_x%H(i1, ne, i2, 2 * component - 1) &
                          + b%top_x%H(i1, ne, i2, 2 * component)
            ENDIF
          ELSE
            IF( cf .EQ. 'E'  .OR. cf .EQ. 'e') THEN
              WRITE(10,*)   b%bas_y%E(i1, ne, i2, 2 * component - 1) &
                          + b%bas_y%E(i1, ne, i2, 2 * component)
            ELSE
              WRITE(10,*)   b%bas_y%H(i1, ne, i2, 2 * component - 1) &
                          + b%bas_y%H(i1, ne, i2, 2 * component)
            ENDIF
          ENDIF
        ENDDO
      ENDDO
    ELSEIF (ne > g%nyyh) THEN
      DO i2 = g%nxl-b%layers, g%nxgh+b%layers, 1
        DO i1 = g%nxl-b%layers, g%nxgh+b%layers, 1
          IF (i1 < g%nxl) THEN
            IF( cf .EQ. 'E'  .OR. cf .EQ. 'e') THEN
              WRITE(10,*)   b%bas_x%E(i1, ne, i2, 2 * component - 1) &
                          + b%bas_x%E(i1, ne, i2, 2 * component)
            ELSE
              WRITE(10,*)   b%bas_x%H(i1, ne, i2, 2 * component - 1) &
                          + b%bas_x%H(i1, ne, i2, 2 * component)
            ENDIF
          ELSEIF (i1 > g%nxyh) THEN
            IF( cf .EQ. 'E'  .OR. cf .EQ. 'e') THEN
              WRITE(10,*)   b%top_x%E(i1, ne, i2, 2 * component - 1) &
                          + b%top_x%E(i1, ne, i2, 2 * component)
            ELSE
              WRITE(10,*)   b%top_x%H(i1, ne, i2, 2 * component - 1) &
                          + b%top_x%H(i1, ne, i2, 2 * component)
            ENDIF
          ELSE
            IF( cf .EQ. 'E'  .OR. cf .EQ. 'e') THEN
              WRITE(10,*)   b%top_y%E(i1, ne, i2, 2 * component - 1) &
                          + b%top_y%E(i1, ne, i2, 2 * component)
            ELSE
              WRITE(10,*)   b%top_y%H(i1, ne, i2, 2 * component - 1) &
                          + b%top_y%H(i1, ne, i2, 2 * component)
            ENDIF
          ENDIF
        ENDDO
      ENDDO
    ELSE
      DO i2 = g%nxl-b%layers, g%nxgh+b%layers, 1
        DO i1 = g%nxl-b%layers, g%nxgh+b%layers, 1
          IF (i1 < g%nxl) THEN      
            IF( cf .EQ. 'E'  .OR. cf .EQ. 'e') THEN
              WRITE(10,*)   b%bas_x%E(i1, ne, i2, 2 * component - 1) &
                          + b%bas_x%E(i1, ne, i2, 2 * component)
            ELSE
              WRITE(10,*)   b%bas_x%H(i1, ne, i2, 2 * component - 1) &
                          + b%bas_x%H(i1, ne, i2, 2 * component)
            ENDIF
          ELSEIF (i1 > g%nxyh) THEN
            IF( cf .EQ. 'E'  .OR. cf .EQ. 'e') THEN
              WRITE(10,*)   b%top_x%E(i1, ne, i2, 2 * component - 1) &
                          + b%top_x%E(i1, ne, i2, 2 * component)
            ELSE
              WRITE(10,*)   b%top_x%H(i1, ne, i2, 2 * component - 1) &
                          + b%top_x%H(i1, ne, i2, 2 * component)
            ENDIF
          ELSEIF (i2 < g%nzl) THEN
            IF( cf .EQ. 'E'  .OR. cf .EQ. 'e') THEN
              WRITE(10,*)   b%bas_z%E(i1, ne, i2, 2 * component - 1) &
                          + b%bas_z%E(i1, ne, i2, 2 * component)
            ELSE
              WRITE(10,*)   b%bas_z%H(i1, ne, i2, 2 * component - 1) &
                          + b%bas_z%H(i1, ne, i2, 2 * component)
            ENDIF
          ELSEIF (i2 > g%nzyh) THEN
            IF( cf .EQ. 'E'  .OR. cf .EQ. 'e') THEN
              WRITE(10,*)   b%top_z%E(i1, ne, i2, 2 * component - 1) &
                          + b%top_z%E(i1, ne, i2, 2 * component)
            ELSE
              WRITE(10,*)   b%top_z%H(i1, ne, i2, 2 * component - 1) &
                          + b%top_z%H(i1, ne, i2, 2 * component)
            ENDIF
          ELSE
            IF( cf .EQ. 'E'  .OR. cf .EQ. 'e') THEN
              WRITE(10,*) g%E(i1, ne, i2, component)
            ELSE
              WRITE(10,*) g%H(i1, ne, i2, component)
            ENDIF
          ENDIF
        ENDDO
      ENDDO
    ENDIF
  ENDIF

  IF ( co .EQ. 'z' ) THEN
    WRITE(10,*) '% xlabel="x"'
    WRITE(10,*) '% ylabel="y"'
    
    DO i2 = g%nxl-b%layers, g%nxgh+b%layers, 1
      DO i1 = g%nxl-b%layers, g%nxgh+b%layers, 1
        IF(i1 < g%nxl) THEN
          IF( cf .EQ. 'E'  .OR. cf .EQ. 'e') THEN
            WRITE(10,*)   b%bas_x%E(i1, i2, ne, 2 * component - 1) &
                        + b%bas_x%E(i1, i2, ne, 2 * component)
          ELSE
            WRITE(10,*)   b%bas_x%H(i1, i2, ne, 2 * component - 1) &
                        + b%bas_x%H(i1, i2, ne, 2 * component)
          ENDIF
        ELSEIF (i1 > g%nxyh) THEN
          IF( cf .EQ. 'E'  .OR. cf .EQ. 'e') THEN
            WRITE(10,*)   b%top_x%E(i1, i2, ne, 2 * component - 1) &
                        + b%top_x%E(i1, i2, ne, 2 * component)
          ELSE
            WRITE(10,*)   b%top_x%H(i1, i2, ne, 2 * component - 1) &
                        + b%top_x%H(i1, i2, ne, 2 * component)
          ENDIF
        ELSEIF (i2 < g%nyl) THEN
          IF( cf .EQ. 'E'  .OR. cf .EQ. 'e') THEN
            WRITE(10,*)   b%bas_y%E(i1, i2, ne, 2 * component - 1) &
                        + b%bas_y%E(i1, i2, ne, 2 * component)
          ELSE
            WRITE(10,*)   b%bas_y%H(i1, i2, ne, 2 * component - 1) &
                        + b%bas_y%H(i1, i2, ne, 2 * component)
          ENDIF
        ELSEIF (i2 > g%nyyh) THEN
          IF( cf .EQ. 'E'  .OR. cf .EQ. 'e') THEN
            WRITE(10,*)   b%top_y%E(i1, i2, ne, 2 * component - 1) &
                        + b%top_y%E(i1, i2, ne, 2 * component)
          ELSE
            WRITE(10,*)   b%top_y%H(i1, i2, ne, 2 * component - 1) &
                        + b%top_y%H(i1, i2, ne, 2 * component)
          ENDIF
        ELSE
          IF (ne < g%nzl) THEN
            IF( cf .EQ. 'E'  .OR. cf .EQ. 'e') THEN
              WRITE(10,*)   b%bas_z%E(i1, i2, ne, 2 * component - 1) &
                          + b%bas_z%E(i1, i2, ne, 2 * component)
            ELSE
              WRITE(10,*)   b%bas_z%H(i1, i2, ne, 2 * component - 1) &
                          + b%bas_z%H(i1, i2, ne, 2 * component)
            ENDIF
          ELSEIF (ne > g%nzyh) THEN
            IF( cf .EQ. 'E'  .OR. cf .EQ. 'e') THEN
              WRITE(10,*)   b%top_z%E(i1, i2, ne, 2 * component - 1) &
                          + b%top_z%E(i1, i2, ne, 2 * component)
            ELSE
              WRITE(10,*)   b%top_z%H(i1, i2, ne, 2 * component - 1) &
                          + b%top_z%H(i1, i2, ne, 2 * component)
            ENDIF
          ELSE
            IF( cf .EQ. 'E'  .OR. cf .EQ. 'e') THEN
              WRITE(10,*) g%E(i1, i2, ne, component)
            ELSE
              WRITE(10,*) g%H(i1, i2, ne, component)
            ENDIF
          ENDIF
        ENDIF
      ENDDO
    ENDDO
  ENDIF

  WRITE(10,*)
  WRITE(10,*) '$ END'
  CLOSE(10)

END SUBROUTINE density_pml

END MODULE debug_pml