Mixed-Level Simulation with CIDER

tuto
2013-01-10
2013-06-12
  • tuto
    tuto
    2013-01-10

    Hi everybody,

    I'm simulating the circuit in the file diode.cir (from the examples in ngspice) using CIDER. I've got currents and voltages without a problem. Now I want to save or plot some internal parameters from the diode, like doping, potential, concentration, etc. This is possible using the output card in CIDER (ngspice manual, page 477). But it's not working. Perhaps something is missing in my script…

    One-Dimensional Diode Simulation

    * Several simulations are performed by this file.
    * They are:
    *   1. An operating point at 0.7v forward bias.
    *   2. An ac analysis at 0.7v forward bias.
    *   3. The forward and reverse bias characteristics from -3v to 2v.

    Vpp 1 0 0.7v (PWL 0ns 3.0v 0.01ns -6.0v) (AC 1v)
    Vnn 2 0 0v
    D1  1 2 M_PN AREA=2

    .model M_PN numd level=1
    + ***************************************
    + *** One-Dimensional Numerical Diode ***
    + ***************************************
    + options defa=1p
    + x.mesh loc=0.0 n=1
    + x.mesh loc=1.3 n=201
    + domain   num=1 material=1
    + material num=1 silicon
    + mobility mat=1 concmod=ct fieldmod=ct
    + doping gauss p.type conc=1e20 x.l=0.0  x.h=0.0 char.l=0.100
    + doping unif  n.type conc=1e16 x.l=0.0  x.h=1.3
    + doping gauss n.type conc=5e19 x.l=1.3  x.h=1.3 char.l=0.100
    + models bgn aval srh auger conctau concmob fieldmob
    + method ac=direct
    + output op.debug doping psi

    .option acct bypass=0 abstol=1e-18 itl2=100
    .dc Vpp -3.0v 2.0001v 50mv
    .control
    save all @d1 @d1
    run diode.raw
    plot i(vnn)
    .endc
    .END

    Thanks in advance

    Tuto

     
  • Goehler
    Goehler
    2013-01-10

    Hi Tuto,

    postprocessing for cider is rather poor. You can plot quantites in a not very nice way if you load a cider solution file into ngspice and type in some plot commands. Best way is to store the solution files and read them out with Scicoslab or Scilab. By means of imagemagick you can also create impressive movies.
    Cider works on regular grids, that simplifies postprocessing dramatically.

    Below is a cir file (2D pin Diode) and a sample Scicoslab script. For movies, you need to install imagemagick.
    The best raw file viewer is SpiceOpus.

    Best regards

    Lutz

    =======================================================

    TWO-DIMENSIONAL PIN-DIODE CIRCUIT

    VIN 100 0 200.0v
    LS1 100 1 1uH IC=10.0
    L1  1   2 500uH IC=10.0
    R1  2   3 20
    VL  3   4 0.0v
    VD  4   5 0.0v
    D1  5   1 M_PIN AREA=66.66 save
    S1  4   0 10 0 SW1 ON
    VS1 10  0 5.0v PULSE(5.0v 0.0v 10ns 1ns 1ns 1us 2us)

    .MODEL M_PIN NUMD LEVEL=2
    + options defw=1000u

    + x.mesh n=1 l=0.0
    + x.mesh n=3 l=100.0
    +
    + y.mesh n=1   l=0.0
    + y.mesh n=3   l=3.0
    + y.mesh n=13  l=5.0
    + y.mesh n=166 l=153.0
    + y.mesh n=176 l=155.0
    + y.mesh n=179 l=158.0
    +
    + domain num=1 material=1 x.l=0.0 x.h=100.0 y.l=0.0 y.h=158.0

    + material num=1 silicon tn=0.5us tp=0.5us
    +
    + electrode num=1 x.l=0.0 x.h=100.0 y.l=0.0 y.h=0.0
    + electrode num=2 x.l=0.0 x.h=100.0 y.l=158.0 y.h=158.0
    +
    + doping gauss p.type lat.unif conc=1.0e20 char.len=1.076 y.h=0.0
    + doping unif  n.type conc=1.0e14
    + doping gauss n.type lat.unif conc=1.0e20 char.len=1.076 y.l=158.0
    +
    + mobility mat=1 concmod=sg fieldmod=sg
    + mobility mat=1 elec major mumax=1000.0 mumin=100.0 ntref=1.0e16 ntexp=0.8 vsat=1.0e7 vwarm=3.0e6
    + mobility mat=1 elec minor mumax=1000.0 mumin=200.0 ntref=1.0e17 ntexp=0.9
    + mobility mat=1 hole major mumax=500.0  mumin=50.0  ntref=1.0e16 ntexp=0.7 vsat=8.0e6 vwarm=1.0e6
    + mobility mat=1 hole minor mumax=500.0  mumin=150.0 ntref=1.0e17 ntexp=0.8
    +
    + models bgn srh auger conctau concmob fieldmob
    +
    + output rootfile="devres"

    .MODEL SW1 sw vt=2.5 vh=0.5 ron=1.0e-6 roff=1.0e6
    .OPTION ACCT BYPASS=1
    .TRAN 10ns 2us UIC

    .END

    =======================================================

    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
    // This script reads does the following:
    // - read a NGSPICE/Cider raw file (Version of NGSpice > 24)
    // - displays Cider results
    // - plots into files
    //
    // (V1.3/06.06.2012/LGO)
    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////

    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
    // Begin: Initialization 1/2
    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////

    // Clear comand window
    clc;

    // Clear all variables
    clear;

    // Clear current graphics
    clf();

    // Close all open graphic windows
    xdel;

    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
    // End: Initialization 1/2
    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////

    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
    // Begin Function: Convert/Read Cider file
    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////

    function  = ReadCiderResultFile(CiderFile,OutMessage)

        // Set number of variables (number of values per row)
        nvars = 0.0;

        // Set number of points
        npoints = 0.0;

        // Set  number of points in x direction
        npointsx = 0.0;

        // Set number of points in y direction
        npointsy = 0.0;

        // Open Cider file
        fd = mopen(CiderFile,'rb');

        // Create header array
        Header = ;

        // Create header line string
        HeaderLine = "";
       
        // Create points matrix
        nM = ;

        // Read header
        HeaderLine = mgetl(fd,1);
        while HeaderLine <> "Binary:",
         Header = ;
         HeaderLine = mgetl(fd,1);
        end
        Header = ;

        // Get number of variables
        StringFoundArray = grep(Header,"No. Variables:");
         = msscanf(Header(StringFoundArray(2)),"%s %s %i");

        // Get number of points
        StringFoundArray = grep(Header,"No. Points:");
         = msscanf(Header(StringFoundArray(2)),"%s %s %i");

        // Get dimensions
        StringFoundArray = grep(Header,"Dimensions:");
         = msscanf(Header(StringFoundArray(1)),"%s %i , %i");
       
        // Create points matrix
        nM = ;

        // Message
        if (OutMessage > 0.5) then
         disp("Cider file: " + CiderFile);
         disp("Variables: " + string(nvars));
         disp("Points(x): " + string(npointsx));
         disp("Points(y): " + string(npointsy));
         disp("Points(total): " + string(npoints));
        end

        // Create result matrix
        RM = ;

        // Create result matrix row
        RMR = ;
       
        // Read in  
        while (meof(fd) == 0)
         // Read all variables for one point
         RMR = mget(nvars,'dl',fd);
         // Update result matrix
         RM = ; 
        // End of while
        end

        // Close Cider file
        mclose(CiderFile);

        // Get dimensions
         = size(RM);

        // Message
        if (OutMessage > 0.5) then
         disp("Done.");
         disp("Result: RM");
         disp("Dimensions: " + string(RMRows) + "x" + string(RMColumns));
         disp("");
        end

    endfunction

    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
    // End Function: Convert Cider file
    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////

    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
    // Begin Function: GetFileNameForTime
    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////

    function  = GetFileNameForTime(CiderLogFile, TargetTime)

        // Set log array
        LogArray = ;

        // Set log array line
        LogArrayLine = ;

        // Open Cider log file
        fdl = mopen(CiderLogFile,'rb');

        // Read in log file
        LogArrayLine = mgetl(fdl,1);
        while (meof(fdl) == 0),
         LogArray = ;
         LogArrayLine = mgetl(fdl,1);
        end
        LogArray = ;
        mclose(fdl);

        // Get log array size
         = size(LogArray);

        // Get time and file name
         = msscanf(LogArray(1),"%s %s %s %s %e");

        // Extract file name
        FileName = "";
        for i = 1:1:(length(FileNameRaw)-1),
         FileName = FileName + part(FileNameRaw,i);
        end

        // Get file name for TargetTime
        // Get file name for time < TargetTime
        TimeArray = ;
        FileNameArray = ;
        i = 0;
        while ((time < TargetTime)&(i + 1 <= LARows)),
         i = i + 1;
          = msscanf(LogArray(i),"%s %s %s %s %e");
         // Extract file name
         FileName = "";
         for j = 1:1:(length(FileNameRaw)-1),
          FileName = FileName + part(FileNameRaw,j);
         end
         TimeArray(1) = time;
         FileNameArray(1) = FileName;
        end

        // Get file name for time >= TargetTime
        if (i < LARows) then
          = msscanf(LogArray(i),"%s %s %s %s %e");
         // Extract file name
         FileName = "";
         for j = 1:1:(length(FileNameRaw)-1),
          FileName = FileName + part(FileNameRaw,j);
         end
         TimeArray(2) = time;
         FileNameArray(2) = FileName;
        end

        // Decide for result time and result file name
        if ((TargetTime-TimeArray(1) < TimeArray(2)-TargetTime)|(TimeArray(2) <= TimeArray(1))) then
         ResultTime = TimeArray(1);
         ResultFileName = FileNameArray(1);
        else
         ResultTime = TimeArray(2);
         ResultFileName = FileNameArray(2);
        end

        // Get file number   
        IndexArrayDot = strindex(ResultFileName,".");
        NumberString = "";
        for i = IndexArrayDot(1)+1:1:IndexArrayDot(2)-1,
         NumberString = NumberString + part(ResultFileName,i)
        end
         = msscanf(NumberString,"%d");

        // Create complete result file name
        CompleteResultFileName = "";
        ArrayBackslash = ;
        ArrayBackslash = strindex(CiderLogFile, "\");
         = size(ArrayBackslash);
        for i = 1:1:ArrayBackslash(ABColumns),
          CompleteResultFileName = CompleteResultFileName + part(CiderLogFile,i);
        end
        CompleteResultFileName = CompleteResultFileName + ResultFileName;
       
        // Get final time
         = msscanf(LogArray(LARows),"%s %s %s %s %e");

    endfunction

    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
    // End Function: GetFileNameForTime
    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////

    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
    // Begin Function: SlicePlotDopNPEyJnJp
    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////

    function Error = SlicePlotDopNPEyJnJp(RM, ResultTime, xcoord, Delta)

        // Calculate factor for lower bound
        fl = 1.0-Delta*100.0;

        // Calculate Factor for upper bound
        fu = 1.0+Delta*100.0;

        // Get dimensions
         = size(RM);
       
        // Set Error
        if ((RMRows >= 1) & (RMColumns >= 1)) then
         // No Error
         Error = 0.0;

         // Define plot vector 1
         PV1 = ;

         // Fill plot vector 1
         for i=1:1:RMRows,
          if ((RM(i,2) > xcoord*fl) & (RM(i,2) < xcoord*fu)) then
           PV1 = [PV1; ];
          end
         end

         if (PV1 <> ) then
          // Plot vector 1
          subplot(2,2,1)
          plot2d(1.0e6*PV1(:,1),PV1(:,2),style=5,logflag="nl");
          GrapficsObject = gce();
          GrapficsObject.children(1).thickness = 2;
          xgrid(3);
          xtitle("Semi-logarithmic plot of doping along y axis at x = " + string(xcoord*1.0e6) + "um (time: "+string(ResultTime)+"s)","y","N");
         else
          disp("Empty plot vector PV1");
         end

         // Define plot vector 2
         PV2 = ;

         // Fill plot vector 2
         for i=1:1:RMRows,
          if ((RM(i,2) > xcoord*fl) & (RM(i,2) < xcoord*fu)) then
           PV2 = [PV2; ];
          end
         end

         if (PV2<>) then
          // Plot vector 2
          subplot(2,2,2)
          plot2d(1.0e6*PV2(:,1),,style=,logflag="nl");
          GrapficsObject = gce();
          GrapficsObject.children(1).thickness = 2;
          GrapficsObject.children(2).thickness = 2;
          xgrid(3);
          xtitle("Semi-logarithmic plot of electron and hole concentration along y axis at x = " + string(xcoord*1.0e6) + "um (time: "+string(ResultTime)+"s)","y","n,p");
         else
          disp("Empty plot vector PV2");
         end

         // Define plot vector 3
         PV3 = ;

         // Fill plot vector 3
         for i=1:1:RMRows,
          if ((RM(i,2) > xcoord*fl) & (RM(i,2) < xcoord*fu)) then
           PV3 = [PV3; ];
          end
         end

         if (PV3 <> ) then
          // Plot vector 3
          subplot(2,2,3)
          plot2d(1.0e6*PV3(:,1),PV3(:,2),style=5,logflag="nn");
          GrapficsObject = gce();
          GrapficsObject.children(1).thickness = 2;
          xgrid(3);
          xtitle("Linear plot of electric field along y axis at x = " + string(xcoord*1.0e6) + "um (time: "+string(ResultTime)+"s)","y","Ey");
         else
          disp("Empty plot vector PV3");
         end

         // Define plot vector 4
         PV4 = ;

         // Fill plot vector 4
         for i=1:1:RMRows,
          if ((RM(i,2) > xcoord*fl) & (RM(i,2) < xcoord*fu)) then
           PV4 = [PV4; ];
          end
         end

         if (PV4 <> ) then
          // Plot vector 4
          subplot(2,2,4)
          plot2d(1.0e6*PV4(:,1),,style=,logflag="nl");
          GrapficsObject = gce();
          GrapficsObject.children(1).thickness = 2;
          GrapficsObject.children(2).thickness = 2;
          xgrid(3);
          xtitle("Semi-logarithmic plot of electron and hole current density along y axis at x = " + string(xcoord*1.0e6) + "um (time: "+string(ResultTime)+"s)","y","jn,jp");
         else
          disp("Empty plot vector PV4");
         end
      
         // Adjust size
         CF = get("current_figure");
         CF.figure_size = ;
        
         // Adjust size for plotting
         set_posfig_dim(CF.figure_size(1),CF.figure_size(2));
        
        else
         // Error
         Error = 1.0;
        end

    endfunction

    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
    // End Function: SlicePlotDopNPEyJnJp
    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////

    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
    // Begin Function: ThreeDPlotDop
    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////

    function Error = ThreeDPlotDop(RM, nM, ResultTime, AutoAdjust, Ebox)

        // Get dimensions
         = size(RM);

        // Set Error
        if ((RMRows >= 1) & (RMColumns >= 1)) then
         // No Error
         Error = 0.0;

         // x grid vector
         xGrid = ;

         // y grid vector
         yGrid = ;

         // z matrix
         zMatrix = ;

         // z matrix row
         zMatrixRow = ;

         // Compose x grid vector
         for i=1:1:nM(3),
          xGrid = ;   
         end

         // Compose y grid vector
         for i=1:1:nM(4),
          yGrid = ;   
         end

         // Compose z vector
         for i = 1:1:nM(3),
          for j = 1:1:nM(4),
           zMatrixRow = ;
          end
          zMatrix = ;
          zMatrixRow = ;
         end

         // 3D plot
         if (AutoAdjust > 0.5) then
          plot3d(xGrid,yGrid,zMatrix,240,60.0,"",);
         else
          plot3d(xGrid,yGrid,zMatrix,240,60.0,"",,Ebox);
         end
         xtitle("3d plot of doping (time: "+string(ResultTime)+"s)","x","y");
        
         // Adjust size
         CF = get("current_figure");
         CF.figure_size = ;
        
         // Adjust size for plotting
         set_posfig_dim(CF.figure_size(1),CF.figure_size(2));
       
        else
         // Error
         Error = 1.0;
        end
     
    endfunction

    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
    // End Function: ThreeDPlotDop
    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////

    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
    // Begin Function: ThreeDPlotPsi
    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////

    function Error = ThreeDPlotPsi(RM, nM, ResultTime, AutoAdjust, Ebox)

        // Get dimensions
         = size(RM);

        // Set Error
        if ((RMRows >= 1) & (RMColumns >= 1)) then
         // No Error
         Error = 0.0;

         // x grid vector
         xGrid = ;

         // y grid vector
         yGrid = ;

         // z matrix
         zMatrix = ;

         // z matrix row
         zMatrixRow = ;

         // Compose x grid vector
         for i=1:1:nM(3),
          xGrid = ;   
         end

         // Compose y grid vector
         for i=1:1:nM(4),
          yGrid = ;   
         end

         // Compose z vector
         for i = 1:1:nM(3),
          for j = 1:1:nM(4),
           zMatrixRow = ;
          end
          zMatrix = ;
          zMatrixRow = ;
         end

         // 3D plot
         if (AutoAdjust > 0.5) then
          plot3d(xGrid,yGrid,zMatrix,240,60.0,"",);
         else
          plot3d(xGrid,yGrid,zMatrix,240,60.0,"",,Ebox);
         end
         xtitle("3d plot of electrostatic potential (time: "+string(ResultTime)+"s)","x","y");
        
         // Adjust size
         CF = get("current_figure");
         CF.figure_size = ;
        
         // Adjust size for plotting
         set_posfig_dim(CF.figure_size(1),CF.figure_size(2));
       
        else
         // Error
         Error = 1.0;
        end
     
    endfunction

    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
    // End Function: ThreeDPlotPsi
    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////

    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
    // Begin Function: ThreeDPlotN
    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////

    function Error = ThreeDPlotN(RM, nM, ResultTime, AutoAdjust, Ebox)

        // Get dimensions
         = size(RM);

        // Set Error
        if ((RMRows >= 1) & (RMColumns >= 1)) then
         // No Error
         Error = 0.0;

         // x grid vector
         xGrid = ;

         // y grid vector
         yGrid = ;

         // z matrix
         zMatrix = ;

         // z matrix row
         zMatrixRow = ;

         // Compose x grid vector
         for i=1:1:nM(3),
          xGrid = ;   
         end

         // Compose y grid vector
         for i=1:1:nM(4),
          yGrid = ;   
         end

         // Compose z vector
         for i = 1:1:nM(3),
          for j = 1:1:nM(4),
           zMatrixRow = ;
          end
          zMatrix = ;
          zMatrixRow = ;
         end

         // 3D plot
         if (AutoAdjust > 0.5) then
          plot3d(xGrid,yGrid,zMatrix,240,60.0,"",);
         else
          plot3d(xGrid,yGrid,zMatrix,240,60.0,"",,Ebox);
         end
         xtitle("3d plot of electron concentration (time: "+string(ResultTime)+"s)","x","y");
        
         // Adjust size
         CF = get("current_figure");
         CF.figure_size = ;
        
         // Adjust size for plotting
         set_posfig_dim(CF.figure_size(1),CF.figure_size(2));
       
        else
         // Error
         Error = 1.0;
        end
     
    endfunction

    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
    // End Function: ThreeDPlotN
    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////

    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
    // Begin Function: ThreeDPlotP
    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////

    function Error = ThreeDPlotP(RM, nM, ResultTime, AutoAdjust, Ebox)

        // Get dimensions
         = size(RM);

        // Set Error
        if ((RMRows >= 1) & (RMColumns >= 1)) then
         // No Error
         Error = 0.0;

         // x grid vector
         xGrid = ;

         // y grid vector
         yGrid = ;

         // z matrix
         zMatrix = ;

         // z matrix row
         zMatrixRow = ;

         // Compose x grid vector
         for i=1:1:nM(3),
          xGrid = ;   
         end

         // Compose y grid vector
         for i=1:1:nM(4),
          yGrid = ;   
         end

         // Compose z vector
         for i = 1:1:nM(3),
          for j = 1:1:nM(4),
           zMatrixRow = ;
          end
          zMatrix = ;
          zMatrixRow = ;
         end

         // 3D plot
         if (AutoAdjust > 0.5) then
          plot3d(xGrid,yGrid,zMatrix,240,60.0,"",);
         else
          plot3d(xGrid,yGrid,zMatrix,240,60.0,"",,Ebox);
         end
         xtitle("3d plot of hole concentration (time: "+string(ResultTime)+"s)","x","y");
        
         // Adjust size
         CF = get("current_figure");
         CF.figure_size = ;
        
         // Adjust size for plotting
         set_posfig_dim(CF.figure_size(1),CF.figure_size(2));
       
        else
         // Error
         Error = 1.0;
        end
     
    endfunction

    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
    // End Function: ThreeDPlotP
    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////

    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
    // Begin Function: ThreeDPlotEMag
    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////

    function Error = ThreeDPlotEMag(RM, nM, ResultTime, AutoAdjust, Ebox)

        // Get dimensions
         = size(RM);

        // Set Error
        if ((RMRows >= 1) & (RMColumns >= 1)) then
         // No Error
         Error = 0.0;

         // x grid vector
         xGrid = ;

         // y grid vector
         yGrid = ;

         // z matrix
         zMatrix = ;

         // z matrix row
         zMatrixRow = ;

         // Compose x grid vector
         for i=1:1:nM(3),
          xGrid = ;   
         end

         // Compose y grid vector
         for i=1:1:nM(4),
          yGrid = ;   
         end

         // Compose z vector
         for i = 1:1:nM(3),
          for j = 1:1:nM(4),
           zMatrixRow = ;
          end
          zMatrix = ;
          zMatrixRow = ;
         end

         // 3D plot
         if (AutoAdjust > 0.5) then
          plot3d(xGrid,yGrid,zMatrix,240.0,80.0,"",);
         else
          plot3d(xGrid,yGrid,zMatrix,240.0,80.0,"",,Ebox);
         end
         xtitle("3d plot of magnitude of electric field (time: "+string(ResultTime)+"s)","x","y");
        
         // Adjust size
         CF = get("current_figure");
         CF.figure_size = ;
        
         // Adjust size for plotting
         set_posfig_dim(CF.figure_size(1),CF.figure_size(2));
       
        else
         // Error
         Error = 1.0;
        end
     
    endfunction

    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
    // End Function: ThreeDPlotEMag
    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////

    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
    // Begin Function: ThreeDPlotJnMag
    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////

    function Error = ThreeDPlotJnMag(RM, nM, ResultTime, AutoAdjust, Ebox)

        // Get dimensions
         = size(RM);

        // Set Error
        if ((RMRows >= 1) & (RMColumns >= 1)) then
         // No Error
         Error = 0.0;

         // x grid vector
         xGrid = ;

         // y grid vector
         yGrid = ;

         // z matrix
         zMatrix = ;

         // z matrix row
         zMatrixRow = ;

         // Compose x grid vector
         for i=1:1:nM(3),
          xGrid = ;   
         end

         // Compose y grid vector
         for i=1:1:nM(4),
          yGrid = ;   
         end

         // Compose z vector
         for i = 1:1:nM(3),
          for j = 1:1:nM(4),
           zMatrixRow = ;
          end
          zMatrix = ;
          zMatrixRow = ;
         end

         // 3D plot
         if (AutoAdjust > 0.5) then
          plot3d(xGrid,yGrid,zMatrix,240.0,80.0,"",);
         else
          plot3d(xGrid,yGrid,zMatrix,240.0,80.0,"",,Ebox);
         end
         xtitle("3d plot of magnitude of electron current density (time: "+string(ResultTime)+"s)","x","y");
        
         // Adjust size
         CF = get("current_figure");
         CF.figure_size = ;
        
         // Adjust size for plotting
         set_posfig_dim(CF.figure_size(1),CF.figure_size(2));
       
        else
         // Error
         Error = 1.0;
        end
     
    endfunction

    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
    // End Function: ThreeDPlotJnMag
    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////

    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
    // Begin Function: ThreeDPlotJpMag
    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////

    function Error = ThreeDPlotJpMag(RM, nM, ResultTime, AutoAdjust, Ebox)

        // Get dimensions
         = size(RM);

        // Set Error
        if ((RMRows >= 1) & (RMColumns >= 1)) then
         // No Error
         Error = 0.0;

         // x grid vector
         xGrid = ;

         // y grid vector
         yGrid = ;

         // z matrix
         zMatrix = ;

         // z matrix row
         zMatrixRow = ;

         // Compose x grid vector
         for i=1:1:nM(3),
          xGrid = ;   
         end

         // Compose y grid vector
         for i=1:1:nM(4),
          yGrid = ;   
         end

         // Compose z vector
         for i = 1:1:nM(3),
          for j = 1:1:nM(4),
           zMatrixRow = ;
          end
          zMatrix = ;
          zMatrixRow = ;
         end

         // 3D plot
         if (AutoAdjust > 0.5) then
          plot3d(xGrid,yGrid,zMatrix,240.0,80.0,"",);
         else
          plot3d(xGrid,yGrid,zMatrix,240.0,80.0,"",,Ebox);
         end
         xtitle("3d plot of magnitude of hole current density (time: "+string(ResultTime)+"s)","x","y");
        
         // Adjust size
         CF = get("current_figure");
         CF.figure_size = ;
        
         // Adjust size for plotting
         set_posfig_dim(CF.figure_size(1),CF.figure_size(2));
       
        else
         // Error
         Error = 1.0;
        end
     
    endfunction

    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
    // End Function: ThreeDPlotJpMag
    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////

    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
    // Begin Function: ContourPlotDop
    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////

    function Error = ContourPlotDop(RM, nM, ResultTime)

        // Get dimensions
         = size(RM);

        // Set Error
        if ((RMRows >= 1) & (RMColumns >= 1)) then
         // No Error
         Error = 0.0;

         // x grid vector
         xGrid = ;

         // y grid vector
         yGrid = ;

         // z matrix
         zMatrix = ;

         // z matrix row
         zMatrixRow = ;

         // Compose x grid vector
         for i=1:1:nM(3),
          xGrid = ;   
         end

         // Compose y grid vector
         for i=1:1:nM(4),
          yGrid = ;   
         end

         // Compose z vector
         for i = 1:1:nM(3),
          for j = 1:1:nM(4),
           zMatrixRow = ;
          end
          zMatrix = ;
          zMatrixRow = ;
         end

         // Contour plot    
         contourf(xGrid,yGrid,zMatrix,12,1:12,"011","",,);
         xtitle("Contour plot of doping (time: "+string(ResultTime)+"s)","x","y");
         xgrid(1);
        
         // Adjust size
         CF = get("current_figure");
         CF.figure_size = ;
        
         // Adjust size for plotting
         set_posfig_dim(CF.figure_size(1),CF.figure_size(2));

        else
         // Error
         Error = 1.0;
        end
     
    endfunction

    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
    // End Function: ContourPlotDop
    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////

    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
    // Begin Function: VectorFieldPlotE
    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////

    function Error = VectorFieldPlotE(RM, nM, ResultTime)

        // Get dimensions
         = size(RM);

        // Set Error
        if ((RMRows >= 1) & (RMColumns >= 1)) then
         // No Error
         Error = 0.0;

         // x grid vector
         xGrid = ;

         // y grid vector
         yGrid = ;

         // fx matrix
         fxMatrix = ;

         // fx matrix row
         fxMatrixRow = ;

         // fy matrix
         fyMatrix = ;

         // fy matrix row
         fyMatrixRow = ;

         // Compose x grid vector
         for i=1:1:nM(3),
          xGrid = ;   
         end

         // Compose y grid vector
         for i=1:1:nM(4),
          yGrid = ;   
         end

         // Compose fx and fy vector
         for i = 1:1:nM(3),
          for j = 1:1:nM(4),
           fxMatrixRow = ;
           fyMatrixRow = ;
          end
          fxMatrix = ;
          fxMatrixRow = ;
          fyMatrix = ;
          fyMatrixRow = ;
         end

         // Vector field plot
         champ1(xGrid,yGrid,fxMatrix,fyMatrix,5.0,,"011");
         xtitle("Vector field plot of electric field (time: "+string(ResultTime)+"s)","x","y");
         xgrid(1);
        
         // Adjust size
         CF = get("current_figure");
         CF.figure_size = ;
        
         // Adjust size for plotting
         set_posfig_dim(CF.figure_size(1),CF.figure_size(2));

        else
         // Error
         Error = 1.0;
        end
     
    endfunction

    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
    // End Function: VectorFieldPlotE
    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////

    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
    // Begin: Initilization 2/2
    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////

    // Create result matrix
    RM = ;

    // Create points matrix
    nM = ;

    // Create error matrix
    Error = ;

    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
    // End: Initilization 2/2
    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////

    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
    // Begin: Main program
    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////

    // Print header
    printf("\n");
    printf(" Cider visualisation (pindiode) (" + date() + ")\n");
    printf(" ===========================================\n");
    printf("\n");

    // Get cider file at a certain time
    = GetFileNameForTime('C:\Users\lgo\Programming\NGSpiceWork\PinDiode\cider.log',1.0);

    // Set plot path
    PlotPath = 'C:\Users\lgo\Programming\NGSpiceWork\PinDiode\';

    // Set plot index
    PlotIndex = 1;

    // Get result matrix read
    = ReadCiderResultFile(CiderFile,1.0);

    // Define x-coordinate
    xcoord = 50.0e-6;

    // Enable DopNPEyJnJp slice plot
    SlicePlotDopNPEyJnJpFlag = 1.0;

    // Enable 3D Dop plot
    ThreeDPlotDopFlag = 0.0;

    // Enable 3D Psi plot
    ThreeDPlotPsiFlag = 1.0;

    // Enable 3D N plot
    ThreeDPlotNFlag = 0.0;

    // Enable 3D P plot
    ThreeDPlotPFlag = 0.0;

    // Enable 3D EMag plot
    ThreeDPlotEMagFlag = 0.0;

    // Enable 3D JnMag plot
    ThreeDPlotJnMagFlag = 0.0;

    // Enable 3D JpMag plot
    ThreeDPlotJpMagFlag = 0.0;

    // Enable contour Dop plot
    ContourPlotDopFlag = 0.0;

    // Enable E vector field plot
    VectorFieldPlotEFlag = 0.0;

    // Enable ps plot to file
    PlotResultsToFilePS = 0.0;

    // Enable gif plot to file
    PlotResultsToFileGIF = 0.0;

    // Enable gif plot to file
    PlotResultsToFileGIFIdx = 0.0;

    // Enable creation of animation
    CreateAnimation = 1.0;

    // Perform DopNPEyJnJp slice plot
    if (SlicePlotDopNPEyJnJpFlag > 0.5) then
    scf(100);
    Error = ;
    if (PlotResultsToFilePS > 0.5) then
      xs2eps(100, PlotPath+'plotSliceDopNPEyJnJp'+string(ResultTime)+'.eps','l');
    end
    if (PlotResultsToFileGIF > 0.5) then
      xs2gif(100, PlotPath+'plotSliceDopNPEyJnJp'+string(ResultTime)+'.gif');
    end
    if (PlotResultsToFileGIFIdx > 0.5) then
      xs2gif(100, PlotPath+'plotSliceDopNPEyJnJp'+string(PlotIndex)+'.gif');
    end
    end

    // Perform 3D Dop plot
    if (ThreeDPlotDopFlag > 0.5) then
    scf(200);
    Error = [Error; ThreeDPlotDop(RM, nM, ResultTime, 1.0, )];
    if (PlotResultsToFilePS > 0.5) then
      xs2eps(200, PlotPath+'plot3DDop'+string(ResultTime)+'.eps','l');
    end
    if (PlotResultsToFileGIF > 0.5) then
      xs2gif(200, PlotPath+'plot3DDop'+string(ResultTime)+'.gif');
    end
    if (PlotResultsToFileGIFIdx > 0.5) then
      xs2gif(200, PlotPath+'plot3DDop'+string(PlotIndex)+'.gif');
    end
    end

    // Perform 3D Psi plot
    if (ThreeDPlotPsiFlag > 0.5) then
    scf(300);
    Error = [Error; ThreeDPlotPsi(RM, nM, ResultTime, 1.0, )];
    if (PlotResultsToFilePS > 0.5) then
      xs2eps(300, PlotPath+'plot3DPsi'+string(ResultTime)+'.eps','l');
    end
    if (PlotResultsToFileGIF > 0.5) then
      xs2gif(300, PlotPath+'plot3DPsi'+string(ResultTime)+'.gif');
    end
    if (PlotResultsToFileGIFIdx > 0.5) then
      xs2gif(300, PlotPath+'plot3DPsi'+string(PlotIndex)+'.gif');
    end
    end

    // Perform 3D N plot
    if (ThreeDPlotNFlag > 0.5) then
    scf(400);
    Error = [Error; ThreeDPlotN(RM, nM, ResultTime, 1.0, )];
    if (PlotResultsToFilePS > 0.5) then
      xs2eps(400, PlotPath+'plot3DN'+string(ResultTime)+'.eps','l');
    end
    if (PlotResultsToFileGIF > 0.5) then
      xs2gif(400, PlotPath+'plot3DN'+string(ResultTime)+'.gif');
    end
    if (PlotResultsToFileGIFIdx > 0.5) then
      xs2gif(400, PlotPath+'plot3DN'+string(PlotIndex)+'.gif');
    end
    end

    // Perform 3D P plot
    if (ThreeDPlotPFlag > 0.5) then
    scf(500);
    Error = [Error; ThreeDPlotP(RM, nM, ResultTime, 1.0, )];
    if (PlotResultsToFilePS > 0.5) then
      xs2eps(500, PlotPath+'plot3DP'+string(ResultTime)+'.eps','l');
    end
    if (PlotResultsToFileGIF > 0.5) then
      xs2gif(500, PlotPath+'plot3DP'+string(ResultTime)+'.gif');
    end
    if (PlotResultsToFileGIFIdx > 0.5) then
      xs2gif(500, PlotPath+'plot3DP'+string(PlotIndex)+'.gif');
    end
    end

    // Perform 3D EMag plot
    if (ThreeDPlotEMagFlag > 0.5) then
    scf(600);
    Error = [Error; ThreeDPlotEMag(RM, nM, ResultTime, 1.0, )];
    if (PlotResultsToFilePS > 0.5) then
      xs2eps(600, PlotPath+'plot3DEMag'+string(ResultTime)+'.eps','l');
    end
    if (PlotResultsToFileGIF > 0.5) then
      xs2gif(600, PlotPath+'plot3DEMag'+string(ResultTime)+'.gif');
    end
    if (PlotResultsToFileGIFIdx > 0.5) then
      xs2gif(600, PlotPath+'plot3DEMag'+string(PlotIndex)+'.gif');
    end
    end

    // Perform 3D JnMag plot
    if (ThreeDPlotJnMagFlag > 0.5) then
    scf(700);
    Error = [Error; ThreeDPlotJnMag(RM, nM, ResultTime, 1.0, )];
    if (PlotResultsToFilePS > 0.5) then
      xs2eps(700, PlotPath+'plot3DJnMag'+string(ResultTime)+'.eps','l');
    end
    if (PlotResultsToFileGIF > 0.5) then
      xs2gif(700, PlotPath+'plot3DJnMag'+string(ResultTime)+'.gif');
    end
    if (PlotResultsToFileGIFIdx > 0.5) then
      xs2gif(700, PlotPath+'plot3DJnMag'+string(PlotIndex)+'.gif');
    end
    end

    // Perform 3D Jp plot
    if (ThreeDPlotJpMagFlag > 0.5) then
    scf(800);
    Error = [Error; ThreeDPlotJpMag(RM, nM, ResultTime, 1.0, )];
    if (PlotResultsToFilePS > 0.5) then
      xs2eps(800, PlotPath+'plot3DJpMag'+string(ResultTime)+'.eps','l');
    end
    if (PlotResultsToFileGIF > 0.5) then
      xs2gif(800, PlotPath+'plot3DJpMag'+string(ResultTime)+'.gif');
    end
    if (PlotResultsToFileGIFIdx > 0.5) then
      xs2gif(800, PlotPath+'plot3DJpMag'+string(PlotIndex)+'.gif');
    end
    end

    // Perform contour Dop plot
    if (ContourPlotDopFlag > 0.5) then
    scf(900);
    Error = ;
    if (PlotResultsToFilePS > 0.5) then
      xs2eps(900, PlotPath+'plotContDop'+string(ResultTime)+'.eps','l');
    end
    if (PlotResultsToFileGIF > 0.5) then
      xs2gif(900, PlotPath+'plotContDop'+string(ResultTime)+'.gif');
    end
    if (PlotResultsToFileGIFIdx > 0.5) then
      xs2gif(900, PlotPath+'plotContDop'+string(PlotIndex)+'.gif');
    end
    end

    // Perform E vector field plot
    if (VectorFieldPlotEFlag > 0.5) then
    scf(1000);
    Error = ;
    if (PlotResultsToFilePS > 0.5) then
      xs2eps(1000, PlotPath+'plotVFE'+string(ResultTime)+'.eps','l');
    end
    if (PlotResultsToFileGIF > 0.5) then
      xs2gif(1000, PlotPath+'plotVFE'+string(ResultTime)+'.gif');
    end
    if (PlotResultsToFileGIFIdx > 0.5) then
      xs2gif(1000, PlotPath+'plotVFE'+string(PlotIndex)+'.gif');
    end
    end

    // Create animation

    // Begin user time definition
    // Define Start time
    StartTime = 1.0e-15;
    // Define End time
    EndTime = FinalTime;
    // Number of frames
    NumberOfFrames = 100;
    // End unser time definition

    // Create graphs
    if (CreateAnimation > 0.5) then
    // Create graphs
    for PlotIndex=1:1:NumberOfFrames,
      // Get file
       = GetFileNameForTime('C:\Users\lgo\Programming\NGSpiceWork\PinDiode\cider.log',StartTime+(PlotIndex-1)*(EndTime-StartTime)/(NumberOfFrames-1));
      // Get result matrix
       = ReadCiderResultFile(CiderFileAnim,0.0);
      // Create figure
      scf(9999);
      // Create plot
      ThreeDPlotN(RMAnim, nMAnim, ResultTimeAnim, 0.0, );
      // Create auxiliary string
      NZeros = int(log10(PlotIndex));
      NZerosMax = int(log10(NumberOfFrames));
      SZeros = '';
      for i=1:1:(NZerosMax-NZeros+1),
       SZeros = SZeros + '0';
      end
      // Save plot to file
      xs2gif(9999, PlotPath+'plot3DPsi'+SZeros+string(PlotIndex)+'.gif');
      // Clear figure
      clf();
    end
    // Clear graph
    xdel;
    end

    // Plot error
    if (norm(Error,2) > 0) then
    disp("Plot error.")
    else
    disp("Plot performed.")
    end

    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
    // End: Main program
    ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////