function dydt = syspde_para_input(n_eq,i_x_v,y1,t,dxi,mx)
% evaluate first and second spatial derivatives for all variables
for k_eq = 1:1:n_eq
   for i_x = i_x_v(k_eq,1):1:i_x_v(k_eq,2)
   	dydx_matrix(k_eq, i_x) = 0.5*( y1(k_eq,i_x+1) - y1(k_eq,i_x-1) )*dxi;
      d2ydx2_matrix(k_eq,i_x) = ( y1(k_eq,i_x+1) - 2.0*y1(k_eq,i_x) + y1(k_eq,i_x-1) )*dxi^2;
   end
end
%  evaluate temporal derivative for all equations at interior nodes
dydt = zeros(n_eq,mx);
%  input heat transfer parameters
		% L = length of rod
		L = 2.0;
		% for copper
		% rho = density [kg/m^3]
      rho = 8933.1;
		% Cp = heat capacity [J/kg/K]
   	Cp = 0.092*1000.0*4.184;
		% k = thermal conductivity [W/m/K]
   	k = 401.0;
   	%  alpha = thermal diffusivity
   	alpha = k/rho/Cp;
      %  heat transfer coefficient in [W/m^2/K]
      h_tran_coeff = 40;
      % diameter in [m]
      diameter = 0.05;
      % Area in [m^2]
      Area = pi*diameter*L;
      % Volume in [m^3]
      Volume = pi/4*diameter^2*L;
      %factor = 2.0*h_tran_coeff*Area/(rho*Cp);
      afac = h_tran_coeff*Area/(rho*Cp*Volume);
      %  Temperature of the surroundings in Kelvin
		Tsurround = 75.0;
      %factor = h_tran_coeff*Area/(rho*Cp*V);
      afac = h_tran_coeff*Area/(rho*Cp*Volume);
      ffac = afac*Tsurround;
% for k_eq = 1
%  perfectly insulated rod
k_eq = 1;
   for i_x = i_x_v(k_eq,1):1:i_x_v(k_eq,2)
      dydt(k_eq,i_x) = alpha*d2ydx2_matrix(k_eq,i_x);
   end
% for k_eq = 2
%  uninsulated rod
k_eq = 2;
   for i_x = i_x_v(k_eq,1):1:i_x_v(k_eq,2)
      dydt(k_eq,i_x) = alpha*d2ydx2_matrix(k_eq,i_x) + ffac - afac*y1(k_eq,i_x);
   end
%  perfectly insulated rod
k_eq = 3;
   for i_x = i_x_v(k_eq,1):1:i_x_v(k_eq,2)
      dydt(k_eq,i_x) = alpha*d2ydx2_matrix(k_eq,i_x);
   end
%  perfectly insulated rod
k_eq = 4;
   for i_x = i_x_v(k_eq,1):1:i_x_v(k_eq,2)
      dydt(k_eq,i_x) = alpha*d2ydx2_matrix(k_eq,i_x);
   end


