Just make the existing code into a procedure with parameters mu, rho, and g. Like this:

restart: Digits:=20:

ThreePlots:= proc(mu, rho, g)
local
     C_0:= 0.4, d:= 0.0002,
     rho_v:= 1000.0, m:= evalf((4/3)*Pi*rho_v*(d/2)^3),   

     F_D:=proc(v::list)                 
     local
          vx:= v[1], vy:= v[2],
          k:= -evalf((3*Pi*mu*d+(Pi/8)*C_0*rho*(d^2)*sqrt(vx^2+vy^2))/m)
     ;      
          [k*vx,k*vy]
     end proc,                                                      

     Eskritt:=proc(vo::list,h)  
     local F:= F_D(vo);
          [vo[1]+h*F[1], vo[2]+h*(F[2]-g)]
     end proc,

     v_0:= 1.0, h:= 0.002,
 
     OnePlot:= proc(phi, N, C)
     local
          v0:=evalf([v_0*cos(phi), v_0*sin(phi)]),
          i, j, v, r
     ;
          v[0]:=v0;
          for i to N do  v[i]:= Eskritt(v[i-1],h)  end do:
          r[0]:=[0,0]:
          for i to N do  r[i]:= r[i-1]+h*v[i-1]  end do:
          plot([[r[j] $j=0..N]], color= C)
      end proc,

      p1:= OnePlot(Pi/4, 65, green),
      p2:= OnePlot(Pi/5, 55, red),
      p3:= OnePlot(Pi/6, 50, blue)
 ;     
     plots:-display(
          [p1,p2,p3], labels=["",""], scaling= constrained
     )
end proc:

ThreePlots(1.8*10^(-5), 1.2, 9.81);

Then call it again with the new values of mu, rho, and g.

You need to specify a value of alpha in the same places as you specified M. I used alpha=1, but the resulting plots were not within your view. You should eliminate the view option until the plot is otherwise as you want.

b:= A[..,3]:
C:= A[..,[1,2,4..11]]:
x:= LinearAlgebra:-LinearSolve(C,b):
x[8];

You need to define the function f like this:

f:= (t,v)-> 0.0207*v^2 - 893.58;

and not like this:

f(t,v):= ....

Perhaps you are confused by the ?piecewise notation in which the integral answer is presented. You can convert it back to abs form, and you can also convert it to ?Heaviside form.

Download convert_abs.mw

The Interp command that I showed for one of your previous questions works just as well in the two-variable case. So here's your problem (iii):

restart:
U:= [0,1,2]:
V:= [x^2+7, x^3+2*x+3, x^3+5]:
Interp(U,V,y) mod 11;

To see the algorithm for doing this, enter

showstat(`mod/Interp/Interp`);

You just need to pay attention to lines 11-26.

G:= GF(7, 3, x^3+x+1):
G:-`/`(G:-ConvertIn(x^2));

Your problem 1 can be solved by a single Maple command:

Interp([0,1,6], [1,5,2], x) mod 7;

It's entered just like your problem from earlier today.

PDE:= diff(r*diff(T(r,t),r),r)/r = diff(T(r,t),t)/beta:
pdsolve(%, build);

Here's a solution that makes use of parallel processors.

PowerSet:= proc(S::set)
local P,x,s;
     P:= {};
     for x in S do  P:= P, Threads:-Seq({s[], x}, s= {P})  end do;
     {P}
end proc;

pdsolve(diff(y(x,t),x$2)+diff(y(x,t),x)/x=diff(y(x,t),t)/alpha);

Here is why your subs did not work. Consider the expression x+2+sin(x+2) as a tree (in the mathematical sense of tree). The root is the whole expression. The first level of branching divides into three nodes: x, 2, and sin(x+2). The x and the 2 are terminal (leaf) nodes. The sin(x+2) splits into two nodes: the sin and the x+2. The sin is a leaf, and the x+2 splits into two leaf nodes. (Pause now to sketch the tree.) The expression which gets substituted for, i.e., the left side of the = in the subs command, has to be a node in that tree. The x+2 only appears as a single node when it is the argument of sin. The other x and 2 get split immediately from the root.

Does that makes sense?

An expression such a 'InverseErf/erf' needs to have both name quotes and unevaluation quotes, with the unevaluation quotes on the outside: '`InverseErf/erf`'(x). Without the name quotes, the expression is interpreted as an unevaluated division. This is not necessarily your only error; it is just the error that I saw immediately.

When you are inside a procedure F, and you want to return 'F'(x), it is preferable to return 'procname'(x).

It is not the same error. You are much closer to the answer now.

The current problem is that points1[2..6] contain a huge number of integrals that cannot be numerically evaluated because M has no value. In the commands eval(mH, {r= R, epsilon= eps1[i]})  you also need to specify a value for M.

Above the line b:= 2*Pi, add a line a:= whatever. Change the 2*Pi to whatever you want. Change every other occurence of b to (b-a).