linear(vds,vgsvth) := 2*vgsvth*vds-vds^2;
saturation(vds,vgsvth) := vgsvth^2;
draincurrent(vds,vgsvth) := if (vds>vgsvth) then saturation(vds,vgsvth) else linear(vds,vgsvth);
limit(vds) := vds^2;
load(draw);
draw2d(terminal=svg,file_name="IvsV_mosfet",dimensions=,
  color=black,font="Times",nticks=75,proportional_axes=none,font_size=12,
  user_preamble="set termoption enhanced;set mxtics; set mytics; set key bottom right",
  xtics_axis=false,ytics_axis=false,grid=true,points_joined=true,yrange=,
  xaxis=false,yaxis=true,axis_top=true,axis_bottom=true,axis_left=true,axis_right=true,
  line_width=2,xlabel="Drain to source voltage ",
  ylabel="Drain current ",
  label(, ),
  label_alignment=right,label(,,
    ,,
    ,,),
  color=blue,key="",explicit(draincurrent(x,1),x,0,10),explicit(draincurrent(x,2),x,0,10),
  explicit(draincurrent(x,3),x,0,10),explicit(draincurrent(x,4),x,0,10),
  explicit(draincurrent(x,5),x,0,10),explicit(draincurrent(x,6),x,0,10),
  explicit(draincurrent(x,7),x,0,10),
  color=red,line_width=1,line_type=dots,nticks=10,explicit(limit(x),x,0,10)
 );