0
$\begingroup$

I am computing the thermal conductance for the structure 'Al2O3-bilayer graphene-Al2O3' along the z axis, with the following input file.

variable x equal 47.59
variable y equal 41.214148966101
variable z equal 136.636
variable t equal 300.0

#Setup parameters
units metal
atom_style atomic
read_data atomic_structure_atomic
mass 1 16.0
mass 2 27.0
mass 3 12.0
mass 4 12.0

pair_style hybrid/overlay lj/cut 12.0 morse 6.287 tersoff shift -0.0462 kolmogorov/crespi/z 20.0
pair_coeff 1 1 lj/cut 0.00844 3.541
pair_coeff 1 2 lj/cut 0.0135951 4.02
pair_coeff 1 3 lj/cut 0.0 1.0
pair_coeff 1 4 lj/cut 0.0 1.0
pair_coeff 2 2 lj/cut 0.0218989 4.499
pair_coeff 2 3 morse 0.4691 1.738 2.246
pair_coeff 2 4 morse 0.4691 1.738 2.246
pair_coeff * * tersoff ./Copt.tersoff NULL NULL C C
pair_coeff 3 4 none
pair_coeff 3 4 kolmogorov/crespi/z /public23/home/a21000018/lammps/lammps-2Aug2023/potentials/CC.KC NULL NULL C C

neighbor 0.3 bin
neigh_modify delay 0 every 1

#energy minimization
fix 5 all box/relax iso 0.0 vmax 0.001
minimize 1.0e-4 1.0e-6 1000 100000
unfix 5

#layers for heat flux
region hot block INF INF INF INF 28.69356 42.12042884
region cold block INF INF INF INF 98.51465116 111.50565116
compute Thot all temp/region hot
compute Tcold all temp/region cold

velocity all create $t 87287

#1st equilibirum run
fix 1 all nvt temp $t $t 0.1
thermo 1000
run 100000
velocity all scale $t
unfix 1
reset_timestep 0

#2nd equilibrium run
fix 1 all nve
fix hot all heat 1 1.0 region hot
fix cold all heat 1 -1.0 region cold

#thermal conductivity calculation
compute ke all ke/atom
variable temp atom c_ke/1.5/8.617333262145e-5
compute layers all chunk/atom bin/1d z lower 0.02 units reduced  #1/0.02 chunks
fix 2 all ave/chunk 10 100 1000 layers v_temp file profile.heat ave running

run 200000

I set the vacuum zone with the thickness about 27 Angstroms on both sides of Al2O3 along the z axis. I set Al2O3 with the thickness of about 13.0 angstroms as the heat and cold bath on each side respectively. Then, I divide the system along the z axis into 50 chunks and compute the average temperature for each chunk along the z axis.

After the calculation, I obtain the following profile.heat file.

# Chunk-averaged data for fix 2 and group all
# Timestep Number-of-chunks Total-count
# Chunk Coord1 Ncount v_temp
1000 50 19600
  1 0.01 0 0
  2 0.03 0 0
  3 0.05 0 0
  4 0.07 0 0
  5 0.09 0 0
  6 0.11 0 0
  7 0.13 0 0
  8 0.15 0 0
  9 0.17 0 0
  10 0.19 276 1.07973
  11 0.21 692.12 1.03563
  12 0.23 693.95 1.35079
  13 0.25 695.36 1.90038
  14 0.27 711.54 2.85858
  15 0.29 707.26 4.93062
  16 0.31 700.61 7.99583
  17 0.33 680.42 12.3902
  18 0.35 699.79 21.705
  19 0.37 677.82 35.7687
  20 0.39 673.49 100.451
  21 0.41 751.61 496.006
  22 0.43 738.44 947.772
  23 0.45 630.25 1782.94
  24 0.47 331.65 2886.43
  25 0.49 146.74 4356.89
  26 0.51 43.7 3704.75
  27 0.53 16.71 487.923
  28 0.55 132.05 1498.6
  29 0.57 468.22 978.695
  30 0.59 737.44 1035.62
  31 0.61 840.58 903.464
  32 0.63 735.25 550.567
  33 0.65 648.32 87.6685
  34 0.67 682.98 35.6915
  35 0.69 687.39 21.289
  36 0.71 699 12.4074
  37 0.73 697.2 7.51011
  38 0.75 693.56 4.6302
  39 0.77 709.93 2.9444
  40 0.79 685.95 1.96727
  41 0.81 702.72 1.44519
  42 0.83 599.27 1.40809
  43 0.85 12.68 1.1996
  44 0.87 0 0
  45 0.89 0 0
  46 0.91 0 0
  47 0.93 0 0
  48 0.95 0 0
  49 0.97 0 0
  50 0.99 0 0
2000 50 39200
  1 0.01 0 0
  2 0.03 0 0
  3 0.05 0 0
  4 0.07 0 0
  5 0.09 0 0
  6 0.11 0 0
  7 0.13 0 0
  8 0.15 0 0
  9 0.17 0 0
  10 0.19 277.695 1.30917
  11 0.21 694.195 1.26037
  12 0.23 695.34 1.57993
  13 0.25 693.59 2.18105
  14 0.27 712.74 3.43832
  15 0.29 705.385 5.59761
  16 0.31 702.29 9.25588
  17 0.33 679.12 14.7859
  18 0.35 696.84 25.0508
  19 0.37 681.95 40.8758
  20 0.39 669.825 117.341
  21 0.41 745.185 567.562
  22 0.43 740.99 1021.9
  23 0.45 646.515 1903.92
  24 0.47 322.755 2986.39
  25 0.49 137.395 4670.12
  26 0.51 49.38 4098.48
  27 0.53 17.1 786.207
  28 0.55 138.88 1544.02
  29 0.57 468.29 1069.66
  30 0.59 743 1153.05
  31 0.61 831.76 979.571
  32 0.63 732.585 580.876
  33 0.65 641.46 105.548
  34 0.67 679.785 43.4096
  35 0.69 687.725 26.02
  36 0.71 696.585 15.4842
  37 0.73 702.685 9.30317
  38 0.75 694.115 5.53134
  39 0.77 711.715 3.48605
  40 0.79 684.45 2.28466
  41 0.81 701.345 1.69112
  42 0.83 603.35 1.69335
  43 0.85 13.98 1.49991
  44 0.87 0 0
  45 0.89 0 0
  46 0.91 0 0
  47 0.93 0 0
  48 0.95 0 0
  49 0.97 0 0
  50 0.99 0 0
......
......
......
200000 50 3920000.0000000005
  1 0.01 433.277 2231.59
  2 0.03 430.207 2250.46
  3 0.05 428.257 2283.93
  4 0.07 429.139 2327.68
  5 0.09 426.715 2382.88
  6 0.11 420.018 2475.78
  7 0.13 415.449 2584.74
  8 0.15 410.974 2694.39
  9 0.17 410.838 2814.73
  10 0.19 421.141 2976.5
  11 0.21 420.503 3140.34
  12 0.23 409.841 3331.33
  13 0.25 399.789 3518.87
  14 0.27 389.668 3741.35
  15 0.29 389.307 3988.3
  16 0.31 393.047 4220.35
  17 0.33 392.457 4453.06
  18 0.35 391.273 4672.64
  19 0.37 383.784 4808.32
  20 0.39 367.696 4897.24
  21 0.41 360.046 4983.01
  22 0.43 355.259 5069.2
  23 0.45 344.827 5079.11
  24 0.47 327.386 5034.37
  25 0.49 316.369 5014.72
  26 0.51 311.02 4912.15
  27 0.53 316.006 4774.22
  28 0.55 327.864 4724.34
  29 0.57 342.447 4661.86
  30 0.59 353.307 4623.02
  31 0.61 356.999 4497.71
  32 0.63 357.796 4358.9
  33 0.65 359.213 4176.5
  34 0.67 364.638 4002.55
  35 0.69 369.593 3807.51
  36 0.71 378.131 3633.59
  37 0.73 392.692 3501.1
  38 0.75 404.955 3357.08
  39 0.77 419.801 3176.8
  40 0.79 423.498 3012.68
  41 0.81 426.129 2860.85
  42 0.83 432.35 2717.52
  43 0.85 420.236 2588.23
  44 0.87 415.642 2465.71
  45 0.89 420.318 2392.11
  46 0.91 423.575 2327.77
  47 0.93 426.362 2271.88
  48 0.95 428.324 2238.13
  49 0.97 429.921 2221.57
  50 0.99 431.918 2217.88

I wrote a post-processing code to take the time-averaged temperature value for each chunk; plot these averaged temperature values v.s. the distance along the z axis in order to obtain the temperature gradient. What I find is that the temperature increases from Al2O3 on the hot side to the bilayer graphene in the middle of the system; and then, decreases from the bilayer graphene to the Al2O3 on the cold side.

This is the temperatgure v.s. distance figure. enter image description here

I am puzzled by these temperature gradient figure. I suppose the temperature for these chunks should keep on decreasing from the hot side to the cold side.

In this simulation, the Al2O3 on each side is about 3 unit-cell thick along the z axis and 5 unit-cell thick along both x and y axes. Is it because the Al2O3 is too thin and the thickness is smaller than the mean free path of Al2O3 or something else?

This is my own Fortran code, which processes the profile.heat file and outputs the distance v.s. temperature data.

PROGRAM SLOPE
IMPLICIT NONE

INTEGER, PARAMETER :: dp = SELECTED_REAL_KIND(15,14)
INTEGER :: i, j, iso, ti, cn, an, ln
CHARACTER (LEN=100) :: str
INTEGER, ALLOCATABLE :: ia(:)
REAL (KIND=dp), ALLOCATABLE :: da(:,:,:), sa(:)

OPEN (UNIT=3, FILE='profile.heat', STATUS='OLD', IOSTAT=iso)
OPEN (UNIT=4, FILE='slope.dat', STATUS='UNKNOWN')

DO i =1, 3, 1
   READ (UNIT=3, FMT=*, IOSTAT=iso)
END DO
READ (UNIT=3, FMT=*, IOSTAT=iso) ti, cn, an
DO i =1, cn, 1
   READ (UNIT=3, FMT=*, IOSTAT=iso)
END DO
ln = cn + 1
DO WHILE (iso == 0)
         READ (UNIT=3, FMT=*, IOSTAT=iso)
         ln = ln + 1
END DO
ln = ln - 1
print *,ln
ln = ln / (cn+1)
print *,ln

REWIND (UNIT=3)
ALLOCATE (ia(cn))
ALLOCATE (da(ln,cn,3))
ALLOCATE (sa(cn))
DO i =1, 3, 1
   READ (UNIT=3, FMT=*, IOSTAT=iso)
END DO
DO i = 1, ln, 1
   READ (UNIT=3, FMT=*, IOSTAT=iso)
   DO j = 1, cn, 1
      READ (UNIT=3, FMT=*, IOSTAT=iso) ia(j), da(i,j,:)
   END DO
END DO

sa = 0.0d0
DO i = 1, ln, 1
   DO j = 1, cn, 1
      sa(j) = sa(j) + da(i,j,3)
   END DO
END DO

sa = sa / DBLE(ln)

DO i = 1, cn, 1
   WRITE (UNIT=4, FMT=*) da(1,i,1), sa(i)
END DO

DEALLOCATE (ia)
DEALLOCATE (da)
DEALLOCATE (sa)

CLOSE (UNIT=3)
CLOSE (UNIT=4)

STOP
END PROGRAM SLOPE

Would anyone please take a look at the input file and tell me what is not set properly?

Thank you very much.

$\endgroup$
3
  • 1
    $\begingroup$ Please add the graphs and the code for plotting that you used. $\endgroup$ Commented Apr 18 at 9:05
  • $\begingroup$ @Vandan Revanur Thank you for the reply. Would you please tell me how to upload figure over the stack exchange? Is there any attachment button for uploading the figure and code? Thank you again. $\endgroup$
    – Kieran
    Commented Apr 18 at 13:20
  • $\begingroup$ @Vandan Revanur I upload the figure and my own Fortan code to process the profile.heat file. Would you please have a look at it and give me some suggestions/comments? Thank you. $\endgroup$
    – Kieran
    Commented Apr 18 at 13:53

0

You must log in to answer this question.

Browse other questions tagged .