fdtd simulation, plotting with gnuplot, writting in perl

# 9月13日 于成都黄龙溪
 1 #!/usr/bin/perl
 2 
 3 # Author : Leon  Email: yangli0534@gmail.com
 4 # fdtd simulation , plotting with gnuplot, writting in perl  
 5 # perl and gnuplot software packages should be installed before running this program
 6 
 7 #use Time::HiRes qw(sleep);
 8 #use autodie qw(:all);
 9 print "@
";
10 my $terminal = "";
11     open GNUPLOT_TERM, "echo 'show terminal;' | gnuplot 2>&1 |";
12     while (<GNUPLOT_TERM>) {
13     if (m/terminal type is (w+)/) {
14         $terminal=$1;
15     }
16     }
17     close GNUPLOT_TERM;
18 
19     # unfortunately, the wxt terminal type does not support positioning. 
20     # hardcode it...
21     $terminal  = "x11";
22 
23 open my $PIPE ,"| gnuplot " || die "Can't initialize gnuplot number 
";
24 
25 print $PIPE "set size 0.85, 0.85
";
26 print $PIPE "set term png size 600, 400
";
27 
28 my $title = "fdtd simulation by leon,yangli0534\\@"."gmail.com";
29 print $PIPE "set terminal gif animate
";# terminal type: png
30 print $PIPE "set output "fdtd_simulation_v0.1.gif"
";#output file name 
31 print $PIPE  "set title "{/Times:Italic $title}"
";# title name and font
32 #print $PIPE  "set title  "fdtd simulation by leon,yangli0534\\@ gmail.com"
";# title name and font
33 print $PIPE  "set title  font ",15"  norotate tc rgb "white"
";
34 print $PIPE "unset key
";
35 print $PIPE "set tics textcolor rgb "white"
";# text color
36 print $PIPE "set border lc rgb "orange"
";
37 print $PIPE "set grid lc rgb"orange"
";
38 print $PIPE "set object 1 rectangle from screen 0,0 to screen 1,1 fc rgb "gray10" behind
";#background color
39 print $PIPE "set xlabel" {/Times:Italic distance: wave length}" tc rgb "white" 
";# xlabel
40 print $PIPE "set ylabel"{/Times:Italic amplitude: v}" tc rgb "white"
";#ylabel
41 print $PIPE "set autoscale
"; 
42 
43 my $size = 400;#physical distance
44 my @ez;#electric field
45 my @hy;#magnetic field
46  
47 my $imp0 = 377.0;
48 #initalization
49 for (my $i = 0; $i < $size; $i++){
50     $ez[$i] = 0;
51     $hy[$i] = 0;
52      
53 } 
54 my $qTime;
55 my $MaxTime = 1850;
56 my $pi = 3.141592563589793;
57 print $PIPE "set xrange [0:$size-1]
";
58 my $mm = 0;
59 
60 #do time stepping
61 for($qTime = 0; $qTime < $MaxTime; $qTime+=5){
62 
63     # update magnetic field
64     for( $mm = 0; $mm < $size - 1; $mm++){
65         $hy[$mm] = $hy[$mm] + ($ez[$mm+1] - $ez[$mm])/$imp0;
66     }
67 
68     # update electric field
69     for( $mm = 1; $mm < $size ; $mm++){
70         $ez[$mm] = $ez[$mm] + ($hy[$mm] - $hy[$mm-1])*$imp0;
71     }
72     
73     if($qTime % 10 == 0){
74 
75         print $PIPE "plot "-" w l  lw 3 lc rgb "green"
";
76         my $cnt = 0;
77         for my $elem ( @ez) {
78         #print " ".$elem;
79             print $PIPE $cnt." ".$elem."
";
80             $cnt += 1;
81         }
82         print $PIPE "e
";
83     } 
84     #hardwire a source
85     $ez[0] = exp(-($qTime - 30.0)*($qTime - 30.0)/100);
86 }
87 
88 #print $PIPE "set terminal x11
";
89 
90 print $PIPE "set output
";
91 
92 close($PIPE);
OPTIMISM, PASSION & HARDWORK