This can be solved with LinearAlgebra:-LinearSolve in exactly the same way as a problem that you asked yesterday.

The R5 must be a typo. The given vectors are clearly in R⁴.

There seems to be some missing connection between your two "Suppose...." statements. The only variables in the first are the us, and there are no us in the second. Is A the matrix whose columns are the us? Is v a purely symbolic vector? If both of these are "yes", then b1 = u1[1]*v[1] + u2[1]*v[2] + u3[1]*v[3] + u4[1]*v[4], or, more succinctly, A[1, ..] . v.

You want the coefficients of a linear combination of vectors that equals another vector. The command for that is LinearSolve. Note that the xs do not play any numeric role in this process. The solution that we obtain will be correct for any xs.

LinearAlgebra:-LinearSolve(
    < <71, 73, -153, -259, -108, 245> |
      <37, 189, 287, -167, 279, -51>  |
      <-199, -200, -62, 59, 262, -70> |
      <48, 295, 18, 235, 209, 279>
    >,
    <6143, 20711, 8974, -30368, 18964, 17937>
);

The commands

A:= <
    < -17, -29, 6, -81, 20>     |
    <-65, -11, -47, 18, -15>    |
    <-240, 90, -265, 495, -175> |
    <-53, 70, 84, -80, 61>      |
    <9, 0, 46, -55, -37>        |
    < 176, -280, -520, 540, -96>
>;
LinearAlgebra:-NullSpace(A);

will return a set of two vectors. Any linear combination of those vectors other than identical 0 is a solution to the problem. The solution can also be obtained from LinearSolve:

LinearAlgebra:-LinearSolve(A, <0, 0, 0, 0, 0>);

The solution will be returned as a single vector with two parameters (with names such as _t1[3]). Giving those parameters any values other than both 0 gives a solution to the problem.

Like this:

Has:= (expr, excp, x)-> local Z; has(subsindets(expr, excp, ()-> Z), x):
expr:=1/exp(z)*arcsinh(x*exp(C[1]))+x*C[1]*sin(exp(x))+3*exp(C[1]*y)*sqrt(sin(exp(3*C[1]))):
expr2:=1/exp(z)*arcsinh(x*exp(C[1]))+x*sin(exp(x))+3*exp(C[1]*y)*sqrt(sin(exp(3*C[1]))):
Has(expr, specfunc(exp), C[1]);
                              true

Has(expr2, specfunc(exp), C[1]);
                             false

It's a mathematical convention that when an associative operator with identity is applied to 0 operands (i.e., an empty list of operands), then the return value is the identity of the operator. Compare the following in Maple:

`+`(), `*`(), `union`(), `and`(), `or`();
             0, 1, { }, true, false

The integral int(Dirac(f(x)), x= a..b) simply counts the roots of f(x) in the interval. Assuming that the roots are simple, f(x) is analytic, and a bound dI can be given such that it's known that the box z= a-dI*I..b+dI*I contains only real roots, then the roots can be counted by this simple numeric integration procedure:

RealRootCount:= proc(
    f::algebraic, X::name=range(realcons), DI::realcons:= 0
)
description `Counts roots of analytic function with multiplicity`;
local 
    x:= lhs(X), a:= lhs(rhs(X)), b:= rhs(rhs(X)),
    F:= unapply(diff(f,x)/f/2/Pi, x),
    dI:= `if`(DI=0, (b-a)/1e3, DI), opts:= ('epsilon'= 1e-3, _rest) 
;
    (round@`simplify/zero`@fnormal@evalf[trunc(evalhf(Digits))])( 
        I*Int(F(x+dI*I) - F(x-dI*I), x= a..b, opts) -
        Int(F(a+x*I) + F(b+x*I), x= -dI..dI, opts)      
    )
end proc
:
RealRootCount(sin(x), x= -3/2*Pi..3/2*Pi);
                               3

One-point assumptions such as variable= value are ineffective whether using assume or assuming. Perhaps a warning should be issued when they are attempted.

Instead, do

sol:= pdsolve(eval(pde, alpha= 0), build)

Your code is massively inefficient. One of many possible replacements is

T1:= table(sparse, L1=~[$1..nops(L1)]):
Threads:-Map[2](index, T1, L2);
[2, 3, 8, 12, 15, 0, 0, 0, 0, 0]

subsindets(expr, specfunc(exp), ()-> Z)

IMO, subsindets and evalindets are the most important commands for symbolic/syntactic processing.

If you want to sort^[*1] matrix A by its column C, then use

A[sort(A[.., C], 'output'= 'permutation')];

It doesn't matter whether the matrix is square.

[*1] To be more formal, I'd say "to sort the rows of matrix A using column C as the key, ...."

Ah, this procedure listpart was the very first thing that caught my attention as needing improvement in the lengthy worksheet that you sent me a few days ago. Assuming that nchunk is the number of chunks and that you'd prefer them as equally sized as possible, here's how I'd do it:

ListPartition:= (L::{list, set}, nchunk::posint)->
  local 
        j::thread_local,
        J:= Iterator:-SplitRanks(nops(L), 'numtasks'= nchunk)
    ;
    [Threads:-Seq](L[j[1]..add(j)-1], j= J)        
:
ListPartition([$1..30], 7);
    [[1, 2, 3, 4, 5], [6, 7, 8, 9, 10], [11, 12, 13, 14], 
      [15, 16, 17, 18], [19, 20, 21, 22], [23, 24, 25, 26], 
      [27, 28, 29, 30]]

Why does this bug occur? Somewhere in some child of radnormal, the large integer shown is factored. Here's its true prime factorization:

n:= 646162507019111437893207695980096110233782566593779:
CodeTools:-Usage(ifactor(n));
memory used=7.02MiB, alloc change=256.00KiB, cpu time=109.00ms, real time=142.00ms, gc time=0ns

(5415624023749) (70334392823809) (168749965921) (10052678938039)

We can see that Maple had no trouble (timewise) doing that. However, ifactor has an option easy (see ?ifactor) for very quick usage which gives just "small" prime factors multiplied by some variables that represent composite numbers. Some quick testing shows that the cutoff for what's considered "small" is about 1.3 million.

So let's see what this option does with our n:

CodeTools:-Usage(ifactor(n, easy));
memory used=24.45KiB, alloc change=0 bytes, cpu time=0ns, real time=2.00ms, gc time=0ns

                         _c27_2 _c25_2

This says that it knows that n is a product of two composite numbers---one with 27 digits and one with 25 digits---and it refuses to do more work than that under option easy. The _2s are just session-dependent serial numbers; you can see the _25 and _27 digit counts in your error message with different serial numbers.

Some child of radnormal must've called ifactors(n, easy) but was unprepared to deal with those symbolic composites.

The code of radnormal and its children is written in Maple (not built-in to the kernel).

An expression e(...) with no space after the e (in 2d Input) is processed as a function. So, in your examples, x is treated as a function symbol. When a number appears as a function symbol, it's processed as a constant function, which is why you get 1.

You were only a few steps away from the solution. Note that the GAMMA function is all uppercase and Pi is a predefined constant.

My style is a little different than Rouben's or Tom's with regard to := vs. =. I tend to prefer using fewer := and instead use eval[recurse] to consolidate everything in the dsolve step. But there's nothing wrong with their styles either!

ODEs:= 
    diff(U1(x),x)=-diff(phi(x),x)/(U1(x)-T/U1(x)),
    diff(phi(x),x$2)=
        (1+A1*phi(x)+A2*phi(x)^(3/2)+A3*phi(x)^2)-(M1/U1(x))
:
Evals:= {
    A1 = (2*k-1)/(2*k-3),
    A2 = 8*sqrt(2/Pi)*(beta-1)*k*GAMMA(k)/
        (3*GAMMA(k-0.5)*(2*k-3)^(3/2)),
    A3 = (4*k^2-1)/(2*(2*k-3)^2),
    M1 = 0.1+sqrt(T+(1/A1))
}:
ICs:= U1(0) = M1, phi(0) = 0, D(phi)(0) = 0.001:
Sol:= dsolve(
    eval[recurse]({ODEs, ICs}, Evals),
    numeric, parameters= [T, beta, k]
):
Sol(parameters= [T= 0.1, beta= 0.6, k= 3.5]);
plots:-odeplot(Sol, [x, phi(x)], x= 0..22);
plots:-odeplot(Sol, [x, U1(x)], x= 0..22);
plots:-odeplot(Sol, [x, diff(phi(x),x)], x= 0..22);