I do not understand why you talk about the Physics package.

Is it this package that you really want to use (seems to me a little be complicated for the question you ask).
If you want something else why not

restart
interface(version)
Standard Worksheet Interface, Maple 2015.2, Mac OS X, December 

   21 2015 Build ID 1097895
with(VectorCalculus):
v := <c, d>;  # the display is smarter in a worksheet
           v := c*e[x]+d*e[y]
position := <p, q > + int(v, s)
           position := (c*s+p)*e[x]+(d*s+q)*e[y]

Without any package;

restart
v := <c, d>;
p := t -> <x(t), y(t)>:
         
position := eval(p(t), dsolve({ diff~(p(t), t)=v, p(0)=<alpha, beta> }))
          

Concerning "_i"

You can change _i by _WhatYouWant and get exactly the same kind of results (keep in mind what I told you yb a prevo-ious question about names beginning with an underscore https://www.mapleprimes.com/questions/235461-What-Is-The-Type--obtained-For)
_i has absolutely nothing particular in general:

about(Pi)
Pi:
  is assumed to be: Pi
  a known property having {irrational} as immediate parents
  and {BottomProp} as immediate children.
  mutually exclusive with {0, 1, integer, prime, composite, fraction, rational}

about(_i)
_i:
  nothing known about this object

When you use the Physics package things may change.
See  < Overview of the Physics[Vectors] Subpackage > :

Regarding projected vectors, the Vectors subpackage is designed to work only with cartesian, cylindrical and spherical orthonormal basis and the related systems of coordinates (see Identify), according to
           (_i, _j, _k)
          ...
and further
NOTE: these variables x, y, z, rho, phi, r and theta, as well as _i, _j, _k, _rho, _phi, _r and _theta, respectively used to represent the coordinates and the unit vectors, are automatically protected when the Physics[Vectors] subpackage is loaded.


Concerning "D"

# here f is explicitely defined as a function a a single argument u
f := u -> u*v:
D(f)
                          u -> v

# nothing but the product rule you learnt at school
# https://en.wikipedia.org/wiki/Product_rule
g := u*v:
D(g)
                        D(u) v + u D(v)

Finally, the correct syntax is not D(f(u)) but

D(f)(u)
                               v
D(f)(t)
                               v

Look to D help page to see how to use D correctly
When tou write D(f(u)) the first step is the evaluation of f(u), so it is  D(u*v) that you compute, and thus the result you got.

Question 1
Given two vectors A and B of same type (both row vectors or both column vectors) and same "length", and given two numbers p and q then

C := p*A + q*B

buids a vector C of the same type and length whose element C[i] = p*A[i]+q*B[i]

A := <1 | 2>:
B := <2 | 1>:
C := 2*A + 3*B;
           C :=  [ 8 7 ]

For other componentwise operations on vectors or matrices use ~:

E := A *~ B
            E := Vector[row](2, {(1) = 2, (2) = 2})
F := cos~(A) *~ exp~(B)
            F := Vector[row](2, {(1) = cos(1)*exp(2), (2) = cos(2)*exp(1)})

Question 2
About _names the things are very clear: never ever use such names.
If you are fan with underscores as prefix, just double them:

__a := 3:
__b := 4:
__c := __a + __b;
                               7

This will spare you a lot of problems.
Even if a single underscore can lead to no error at all, your are never sure of that:

# all is fine here
_toto := 3:
_toto^2
                               9

Question 3
Look to the type help page and click on +:

• type/`+` - check for an expression of type `+`
• type/`*` - check for an expression of type `*`
• type/`^` - check for an expression of type `^`

The fact c:=a+b is of  `+` type comes from the way an algebraic operation (here a sum) is represented by a syntactic graph:

c := a+b:
whattype(c)
                               +
ToInert(c)
         _Inert_SUM(_Inert_NAME("a"), _Inert_NAME("b"))

... unless you have a problem with your version.
(BTW your worksheet doesn't contain any plot command: wher do your plots come from ???)

In the RootVector help page search for the string plot and do what it's said:

Download PlotRootedVector_mmcdara.mw

I guess the ambuity comes from the use of the names "create" and "indexed".

Let's take the case of a procedure through two examples. 
In the first one procedure f assigns to a variable a the value 1, but Maple returns a warning telling you it considers a as a local variable. In the second example a is explicitely declared local

f := proc()
a := 1:
end proc
Warning, `a` is implicitly declared local to procedure `f`


interface(warnlevel=0): 
a:=2: 
f := proc() 
a := 1: 
end proc: 
f(); 
a;
                               1
                               2

f := proc()
  local a:
  a := 1:
end proc

a := 'a'; # 
f(); 
# a is still unknown out of procedure f
a; 
# a way to assign 
assign(a=f()): 
a;
                               a
                               1
                               a
# and another one
print():
a:='a'; 
a:=f(): 
a;              
                               a
                               1
                               a

From help page < local > ... If present, specifies the local variables reserved for use by the procedure...
From help page < assign >... assign(a, B) and assign(a = B) commands make the assignment a := B and return NULL.

assign associates a name (a) with the value of an expression (B).

About the lines a:='a' the mechanism of this clearing operation is illustrated here

restart:
anames(user);
a := 1;
anames(user);
                               1
                               a
a := 'a':
anames(user);  # nothing returned, a no longer exists

The use of local (or global too)  enable declaring names, not really (IMO) creating variables ; a name can be associated with an expression and it's common to say this name is a variable
(
I don't know is this reference can be found freely on the web "Introduction to Maple, André Heck, Springer-Verlag, 1993" ?
If it is so, from page 63: "With the procedure assign you can get the status of individual variables".
Then assign, or :=, both give to a name the status of variable ; I guess one can say that assign creates variables.
)

According to the book "Maple V, language Reference Manual, Char et al., 1991", the term variable is used quite lately in it, mainly in refereces to the previous clearing process and he use of the global statement.

global acts basically as local does, excepted that the use of global a; makes a to be known at all levels in your worksheet:

restart:

f := proc()
global a:
a:=1:
end proc:

anames(user);
f():
a;
                               f
                               1

From verion 2019 (I think) it is no monger possible to use the statement global at the upper level of the worksheet, that is to write something like this

restart:
global a:
                               f
                               1

So let me rewrite your question as  "Is it at all possible to assign indexed names"
The answer is yes. From < indexed > help page:

• An indexed expression of the form name[expression sequence] represents a selection operation.
• An indexed name is a valid name.  In particular, since A[1,2,3] is a valid name, so are A[1,2,3][x,y] and A[1,2,3][x,y][2,1].
• The use of the indexed name b[1] does not imply that b is an array.  For instance, if b is an unassigned name, a := b[1] + b[2] + b[1000] simply forms a sum of three indexed names.
• If an assignment is made to name, then name[expression sequence] represents a selection from name using the index sequence expression sequence.  For more information, see selection.

Here are a few valid names:

• R1, R2, R3, ...
• R[1], R[2], R[3], ...

Only the last ones are indexed names:

restart:
R[1] := 1:
R1   := 2:
R[1], R1;  # R[1] is not R1 !!!
print("--------------------------"):

restart:
type(R[1], indexed);
type(R1, indexed);
                              1, 2
                  "--------------------------"
                              true
                             false

An important thing (see bove) is that "The use of the indexed name b[1] does not imply that b is an array" (nor a vector).
The index is not necessarilly an integer:

restart:
type(R[p], indexed);
# but...
type(R[p], indexed[integer]);
                              true
                             false

Be carefull, some names which appear as indexed names are nor indexed names:

restart:
L := [R[1], R1, R__1];
lprint(L);          # look to < lprint > help page
type~(L, indexed);  # look to < ~ > (operators, elementwise (~) ) help page
                        [R[1], R1, R__1]
[R[1], R1, R__1]
                      [true, false, false]

How to assign indexed and non indexed names?
What won't work:

seq(R[n], n=1..3);
f := proc()
  local seq(R[n], n=1..3)
end proc:
                        R[1], R[2], R[3]
Error, `(` unexpected

f := proc()
  local R[1], R[2], R[3]:
end proc;

Error, `[` unexpected

seq(R||n, n=1..3);
f := proc()
  local seq(R||n, n=1..3)
end proc:
                           R1, R2, R3
Error, `(` unexpected

seq(R__||n, n=1..3);
f := proc()
  local seq(R__||n, n=1..3)
end proc:
                        R__1, R__2, R__3
Error, `(` unexpected

And what works

f := proc()
  local R1, R2, R3:
end proc;

f := proc()
  local R__1, R__2, R__3:
end proc;

but neither R1 nor R__1 are indexed names.

If you know the number n of indexed names you want to create 

restart
f := proc(n)
  local R := Vector(n, symbol=R):
  local i:
  for i from 1 to n do R[i]:=i end do:   return R; end proc:

res := f(3); 
            Vector(3, {(1) = R[1], (2) = R[2], (3) = R[3]})

assign(seq(R[i]=res[i], i=1..3)):  # multiple assignment

# check for the values of R[1], R[2] and R[3]
'R[1]' = R[1];
'R[2]' = R[2];
'R[3]' = R[3];
                            R[1] = 1
                            R[2] = 2
                            R[3] = 3

Finally, and it is how I understant you question, you can associate indexed names with expressions at the upper level while doing this

restart:
assign(seq(R[i]=t^(i-1)*cos((i-1)*Pi/2), i=1..3)):

'R[1]' = R[1];
'R[2]' = R[2];
'R[3]' = R[3];
                            R[1] = 1
                            R[2] = 0
                                    2
                           R[3] = -t 

and, more specifically:

restart:

MyFiles := FileTools:-ListDirectory(cat(currentdir()))[1..3];
        MyFiles := ["EgyptianFraction.mw", "matrix_inverse_mmcdara.mw", "Formal_ODE_Solution.mw"]

assign(seq(R[i]=MyFiles[i], i=1..3)):

'R[1]' = R[1];
'R[2]' = R[2];
'R[3]' = R[3];

                  R[1] = "EgyptianFraction.mw"
               R[2] = "matrix_inverse_mmcdara.mw"
                R[3] = "Formal_ODE_Solution.mw"

Since I am not familiar enough with the Physics package, here is an workaround, for what it is worth.

e := expand(lhs(Newton_generalized_MatForm)):
select(has, [op(convert(e, D))], :-D(`ϰ_`)(t)):
convert(coeff(eval(add(%), :-D(`ϰ_`)(t)=J), J), diff)

Nevertheless I can't advise you enough to await a better solution ( @ecterrab ? )

If you want to store several Maple-generated matrices converted into Python code, I propose you this:

restart:
# just a notional example 
with(LinearAlgebra):
t := table([seq(i=RandomMatrix(4, 4), i=1..100)]):

# redirect the default output (terminal) to a file named "python.txt"
f := cat(currentdir(),"/Desktop/python.txt"): 
writeto(f):
for i in sort([indices(t, nolist)]) do
  CodeGeneration:-Python(t[i])
end do:

# reset the default output for subsequent work, if any
writeto(terminal);

python.txt

If you want to read from Maple the Python matrices the file above contains, do this

restart
f    := cat(currentdir(),"/Desktop/python.txt"):
T    := table([]):
line := readline(f):
i    := 1:
while line <> 0 do
   if line <> "> writeto(terminal);" then
     T[i] := convert(parse(StringTools:-Split(StringTools:-Split(line, "(")[-1], ")")[1]), Matrix):
   end if:
   line := readline(f):
   i    := i+1:
end do:
eval(T)

PythonMatrices_mmcdara.mw

It is always possible to rename the vertices as you want:

Download ExtensionField_mmcdara.mw

Generating automatically the vertex names is more complex while not impossible.
Note that non-breaking spaces can also be added in the vertex names for a better rendering and that their colors can be adjusted the same way.

Executing showstat(GroupTheory:-DrawSubgroupLattice) reveals DrawSubgroupLattice offers more capabilities then the corresponding help page (Maple 2015.2) describes (look to SubgroupLattice for more help).

  {
    assignlattice::name := undefined, 
    center := undefined, 
    derived := false, 
    directed::truefalse := false, 
    [edgecolor, edgecolour] := "grey", 
    highlight::{object, set, list({object, set, {object, set, list({object, set})} = anything}), {object, set, list({object, set})} = anything} := [], 
    indices::truefalse := false, 
    labels::identical(:-integers,:-none,:-zuppos,:-ids) := :-integers, 
    normal := true, 
    output::{list(identical(:-graph,:-plot)), 
    identical(:-graph,:-plot)} := :-plot, 
    title := undefined, 
    titlefont := undefined, 
    [vertexcolor, vertexcolour] := "white"
  }

In particular, if you want more flexibility in customizing the graph:

with(GraphTheory):
Gr := DrawSubgroupLattice(g,  'output'=':-graph');
# Then you can use all the features GraphTheory proposes,
# for instance:
Edges(Gr);
GraphTheory:-Vertices(Gr);

Other way:

L  := SubgroupLattice(g);
Gr := convert(L, 'graph')

Even if I guess you have already thought about it

map(D[1], z(u,v)):
type(%, Vector);
      true

diff_mmcdara.mw

For fun (maybe...)

Your polyhedron is a half truncated octahedron 
https://en.wikipedia.org/wiki/Truncated_octahedron
https://mathcurve.com/polyedres/octaedre_tronque/octaedre_tronque.shtml  (in french)

The title refers to the equation its boundary verifies (last reference above)

HalfTruncatedOctahedron.mw

This worksheet also contains a few lines to display all the poyhedra that plots:-polyhedraplot knows.
Could be useful

UPDATED

If it is so, here is way to get it and form a finite difference approximation of the equation 
diff(f(x), x$2)=S(x) when Dirichlet BC hold.

The first part of the worksheet is aimed to derive centered divided difference approximations of Diff(f(x), x) and Diff(f(x), x$2)  (kind of abstract version of fdiff).
The second part of the worksheet is an illutrated aookication of the last result to a numerical approximation of the1D Poisson equation f"(x)-S=0 over (1, 10) with Dirichlet BCs at x=1 and x=10.

The first instance of matrix M is the second-difference matrix of size N by N you are looking for. 
 

Download DD_2.mw

before asking again the same question.
As I said in while answering your previous question, you probably wrote something wrong in the expression you solve for Theta[nu+2] as it refers to values of Theta[n < 0] that are not already defined or to values Theta[1-m] with undefined m.
So this is not even a correct recurrence relation you have to solve.

My feeling is that you use a time-marching scheme (time beung "nu") to solve a partial differential equation whose unknown is Theta(x, nu), starting from the IC Theta(x, 0)=delta(x) and with Neuman or Fourier BCs of convective or radiative type.
If I am right saying this, there are probably better things to do this.


HERE IS YOUR FILE CORRECTED, READ IT... PLEASE
Download ReadMeCarefully.mw

As 2D mode is known as to ease typo errors, here if the 1D mode version

 

Download ReadMeCarefully__1Dmode.mw

Didn't you click on the stop button while the computation of HFS_data_2 and SFS_data_2 were being processed?

I fasten these computations by removing the absolute values (because kappa > 0)
It is likely that further improvements might be made by choosing carefully the integration options.
For instance I used method=_d01akc to compute SFS_data_2 and method=_d01ajc to compute SFS_data_3.
I think it's worth spending some time to try and find the correct method and options.

Here it is: Negativity_(v12_beta_gamma_mu)_mmcdara.mw

(BTW these results have been obtained with Maple 2015.2)

... some corrected, but you must fix yourself a few points:

• Gamma is a symbolic quantity without no given value (not a big issue but the expressions of Theta[nu] remains symbolic and then quite lengthy).
• You  define Theta[nu] (nu =0..8) as a function of Theta[nu-m] without defining m. Same observation than above, but Theta cannot be fully given in explicit form (recurence equation).

The least you must do is to fix m or to define m as a function of nu.

At last, as Theta seems to obey a differential equation, why don't you just solve it?

Download Thermal_Stress_Distribution_mmcdara.mw

To solve the ode I ARBITRARILY considered this was a problem with initial conditions (IC).
If I'm wrong correct this.
Likewise I considered very simple IC, even if your thermal problem suggests that you have convection and radiative conditions.
Here again correct this yourself if you are interesting in dsolving the ode.

I do not believe that Maple can handle the level of abstraction required to derive the results you present.

Nevertheless you could try using the definitions of JacobiP(n, alpha, alpha+1, x) and JacobiP(n, alpha, alpha+2, x)  and, starting from the binomial coefficients these expressions contain, proceed to simplifications in order to express the coefficients of JacobiP(n, alpha, alpha+2, x)  in terms of those of JacobiP(n, alpha, alpha+1, x) :

Download JacobiP.mw