Advanced Search

tomleslie

13876 Reputation

20 Badges

15 years, 181 days
Contact tomleslie

MaplePrimes Activity


These are replies submitted by tomleslie

If your system has a singularity...

tomleslie 13876
March 14 2019
0 0

then it has a singularity. After all, something as simple as 1/(x-1) has an obvious singularity. There is no way to "remove" it

It is *possible* that you might have a floating-point "catastrophic cancellation" problem. You can check this by increasing the precision with which Maple performs floating point calculations (By default, this is 10). Try setting

Digits:=100:

at the top of your worksheet. If you still get a the same error message, then it is (probably) not a floating-point-evaluation issue, but is actually a real singularity - in which case you  are pretty much stuffed.

If the equation system represents some real-world, physical, system, then the occurrence of singularities is rather unusual (but not impossible). Hence if it is a "real-world" system I would suspect a typo somewhere in the definition of the equations you are trying to solve

Have you ever tried reading Maple help?...

tomleslie 13876
March 14 2019
0 0

@izhammulya 
It is very informative and you can actually learn how to solve trivial problems. Since such a simple exercise is beyond you then you could consider the following. I shoul point out that they have very different consequences depending on what you are subsequently going to do with the information provided

> 

restart;
Sol[1] := [72.011122301, [t[4] = -.500000000000000, t[13] = .500000000000000, x[1, 4] = 1, x[1, 13] = 0, x[4, 1] = 0, x[4, 13] = 1, x[13, 1] = 1, x[13, 4] = 0]]:
#
# first way to get any variable from the above
#
  getSol:=z-> select( i->lhs(i)=z, Sol[1][2]):
  getSol(x[1,4]);
  getSol(x[1,13]);
  getSol(x[4,1]);
  getSol(x[4,13]);
  getSol(x[13,4]);

[x[1, 4] = 1]

  

[x[1, 13] = 0]

  

[x[4, 1] = 0]

  

[x[4, 13] = 1]

  

[x[13, 4] = 0]

 (1)
> 

#
# Second way way to get any variable value from the above
#
  map(assign, Sol[1][2]):
  x[1,4];
  x[1,13];
  x[4,1];
  x[4,13];
  x[13,4]

1

  

0

  

0

  

1

  

0

 (2)
> 

 

> 

 


 

Download getVals.mw

The code I supplied...

tomleslie 13876
March 14 2019
0 0

@torabi 

was to solve the specific problem defined in the worksheet you supplied and the pdf you supplied It achieved this

A drawback of asking for a solution to a specific problem is that it is pretty much guaranteed not to work for any other problem. Since you now seem to have a new problem, it is up to you to accurately specify the new problem.

Applying random pieces of code to problems they wer never designed to solve is pretty much guaranteed to fail. Is this difficult to understand?

Because...

tomleslie 13876
March 14 2019
0 0

@a_simsim 

the way you are filling the 3-D array, each entry is essentially a 'pointer' to the same record. In order to get a different record for each array entry, you need to change the command

node[ j , i , k ] := nodeprop :

to

node[ j , i , k ] := copy(nodeprop) :

as in the attached, where the revise code *seems* to give the answers you desire


 

 

 

> 

restart

 

Specify node properties.

 

 

> 

nodeprop := Record(

                    'inode'      = 0 ,
                    'nnb'        = 6 ,
                    'nind'       = Vector( 1 .. 6 ,  0  ) ,
                    'U'          = 100. ,
                    'Ul'         = 100. ,
                    'Uprev'      = 100. ,
                    'Ueq'        = 100. ,
                    'Ueql'       = 100. ,
                    'blin'       = 0.   ,
                    'q'          = 0.   ,
                    'rx'         = 2.   ,
                    'ry'         = 2.   ,
                    'rz'         = 2.   ,
                    'rfin'       = Vector( 1 .. 6 , -1. ) ,
                    'conduct'    = Vector( 1 .. 6 , -1. ) ,
                    'c'          = 2.   ,
                    'tau'        = 0.   ,
                    'k'          = 0.   ,
                    'sumconduct' = 0.    

                  ) :

 

Specify grid dimensions.

 

> 

nx  := 10  :    # X direction dimension.
ny  :=  4  :    # Y direction dimension.
nz  :=  1  :    # Z direction dimension.

 

Establish node array.

 

> 

node := Array( 1 .. ny , 1 .. nx , 1 .. nz ) :

 

Fill up node array with the nodeprop records.

 

> 

for k from 1 to nz do
    for j from 1 to ny do
        for i from 1 to nx do

            node[ j , i , k ] := copy(nodeprop) :

        od : # for i from 1 to nx do
     od : # for j from 1 to ny do
od : # for k from 1 to nz do

 

Fill up inode with values.

 

 

> 

ic := 0 :
for k from 1 to nz do
    for j from 1 to ny do
        for i from 1 to nx do

            ic := ic + 1 :
            
            node[  j , i , k ]:-inode := ic :
            print( j , i , k , node[ j , i , k ]:-inode ) ;

        od : # for i from 1 to nx do
    od : # for j from 1 to ny do
od : # for k from 1 to nz do

1, 1, 1, 1

  

1, 2, 1, 2

  

1, 3, 1, 3

  

1, 4, 1, 4

  

1, 5, 1, 5

  

1, 6, 1, 6

  

1, 7, 1, 7

  

1, 8, 1, 8

  

1, 9, 1, 9

  

1, 10, 1, 10

  

2, 1, 1, 11

  

2, 2, 1, 12

  

2, 3, 1, 13

  

2, 4, 1, 14

  

2, 5, 1, 15

  

2, 6, 1, 16

  

2, 7, 1, 17

  

2, 8, 1, 18

  

2, 9, 1, 19

  

2, 10, 1, 20

  

3, 1, 1, 21

  

3, 2, 1, 22

  

3, 3, 1, 23

  

3, 4, 1, 24

  

3, 5, 1, 25

  

3, 6, 1, 26

  

3, 7, 1, 27

  

3, 8, 1, 28

  

3, 9, 1, 29

  

3, 10, 1, 30

  

4, 1, 1, 31

  

4, 2, 1, 32

  

4, 3, 1, 33

  

4, 4, 1, 34

  

4, 5, 1, 35

  

4, 6, 1, 36

  

4, 7, 1, 37

  

4, 8, 1, 38

  

4, 9, 1, 39

  

4, 10, 1, 40

 (1)

 

> 

node[ 1 , 1 , 1 ]:-inode

1

 (2)

 

Should be 1. Why not ??????

 

 

> 

node[ 2 , 2 , 1 ]:-inode

12

 (3)

 

Should be 12. Why not ??????

 

 

``


 

Download pointer.mw

 

Still don't understand...

tomleslie 13876
March 14 2019
0 0

@izhammulya 

what you are trying to achieve

The attached will take an arbitrary table eg Sol[1], Sol[2],,,,Sol[n] and return the contents by variable name

> 

  restart;
  Sol[1] := [ 72.011122301,
              [ t[4] = -.500000000000000, t[13] = .500000000000000,
                x[1, 4] = 1, x[1, 13] = 0, x[4, 1] = 0, x[4, 13] = 1, x[13, 1] = 1, x[13, 4] = 0
              ]
            ]:
  Sol[2] := [ 53.128387340,
              [ t[6] = -2.00000000000000, t[7] = 0., t[8] = -.999999999999999,
                x[1, 6] = 1, x[1, 7] = 0, x[1, 8] = 0, x[6, 1] = 0, x[6, 7] = 0,
                x[6, 8] = 1, x[7, 1] = 1, x[7, 6] = 0, x[7, 8] = 0, x[8, 1] = 0,
                x[8, 6] = 0, x[8, 7] = 1
              ]
            ]:
  Sol[3] := [ 12.3456,
              [ t[9] = -2.00000000000000, t[10] = 0., t[8] = -.999999999999999,
                z[1, 6] = 1, z[1, 7] = 0, z[1, 8] = 0
              ]
            ]:
  Sol[4] := [ 78.90,
              [ p[9] = -2.00000000000000, p[10] = 0., p[8] = -.999999999999999,
                q[1, 6] = 1, q[1, 7] = 0, q[1, 8] = 0
              ]
            ]:
#
# Procedure which extracts variables by name
# from the supplied table
#
  extractN:= proc( A::table )
                   local varN, getSols;
                   getSols:= (z, zz)-> [ select
                                         ( i-> has(i,z),
                                           { seq
                                             ( zz[j][2][],
                                               j=1..numelems(zz)
                                             )
                                           }
                                         )[]
                                       ];
                   varN:= { seq
                            ( seq
                              ( op( [j,0],
                                    indets( A[k][2])
                                  ),
                                j=1..numelems( A[k][2] )
                              ),
                              k=1..numelems(A)
                            )
                          };
                   return [ seq( getSols(j, A),
                                 j in varN
                               )
                          ]
             end proc:
#
# Organise the content of the supplied table
# by variable name
#
  extractN(Sol);

[[p[8] = -.999999999999999, p[9] = -2.00000000000000, p[10] = 0.], [q[1, 6] = 1, q[1, 7] = 0, q[1, 8] = 0], [t[4] = -.500000000000000, t[6] = -2.00000000000000, t[7] = 0., t[8] = -.999999999999999, t[9] = -2.00000000000000, t[10] = 0., t[13] = .500000000000000], [x[1, 4] = 1, x[1, 6] = 1, x[1, 7] = 0, x[1, 8] = 0, x[1, 13] = 0, x[4, 1] = 0, x[4, 13] = 1, x[6, 1] = 0, x[6, 7] = 0, x[6, 8] = 1, x[7, 1] = 1, x[7, 6] = 0, x[7, 8] = 0, x[8, 1] = 0, x[8, 6] = 0, x[8, 7] = 1, x[13, 1] = 1, x[13, 4] = 0], [z[1, 6] = 1, z[1, 7] = 0, z[1, 8] = 0]]

 (1)
> 

 

> 

 


 

Download sortStuff2.mw

 

The difference...

tomleslie 13876
March 14 2019
0 0

@Hammur 

between x[i] and x__i does not just apply to sums, it is a general notation isssue across the whole of Maple.

x[i] refers to the i-th entry of an indexable quantity eg a list, vector, 1D Array, etc etc

x__i is a simple name (just like xi) exxcep that the 'i' will appear as a subscript in Maple output. It doesn't "index" anything

Observations...

tomleslie 13876
March 13 2019
0 0

The OP's problem can be summarised as:

Given an integer Z, is it possible to construct a right-angled triangle of area Z, all of whose side lengths are rational numbers

This problem is described at

https://en.wikipedia.org/wiki/Congruent_number#Congruent_number_problem

and the (starting) list of integers Z for which it *can* be solved is given at

https://oeis.org/A003273

So far as I can tell, the OP's code functions correctly, in the sense that it returns rational side lengths, which produce the correct area.

His/her problem is execution time.

For the case where the triangle area is 23, execution time for OP's code on my machine is ~6.5 minutes, and no solution is obtained. Note that (according to https://oeis.org/A003273) there should be a solution for Z=23. I can only assume that the required rational number side lengths for this area are really "ugly", and outside the OP's search space.

My conclusions

  1. As wriiten,for the case of Z=23,  the "inner loop" in the OP's code will execute 32,004,000 times, and the if condition it contains will be satisfied 19,455,781 times. That's a lot of execution!
  2. Simply extending the search space for the supplied algorithm in the hope of obtaining a solution will see the execution time increase quadratically, so I'm pretty sure this is not a viable option.
  3. The problem is "parallelizable" so using the Grid toolbox will presumably reduce execution time by a factor equal to the number of threads on your machine Going parallel will result in a linear reduction in execution time while extending the search space will result in a quadratic increase, so this approach would seem to be of little benefit
  4. I reduced the execution time of the OP's algorithm by about 10% with some simple restructuring: not enough to be interesting
  5. Unless someone can come up with a drastically better search algorithm, OP is on a loser :-(

 

I can assure you...

tomleslie 13876
March 12 2019
0 0

This circle will be circular in a Maple worksheet - even if it doesn't look it on this site!

An observation...

tomleslie 13876
March 12 2019
0 0

There is a fundamental difference between a "graph" and a "network".

A "network" must have (at least one) "source" (aka input) and (at least one) "sink" (aka output). A "graph" needs neither.

A "source" can be thought of as a graph vertex all of whose attached edges are directed outwards. A "sink" is a graph vertex all of whose attached edges are directed inwards. This means that a "network"

  1. must be a directed graph,
  2. which in turn means that its associated adjacency matrix cannot be symmetric

You can asses whether an adjacency matrix represents a network rather than a graph by applying the following criteria

  1. A "source" vertex (n) will be characterised by the n-th row in the adjacency matrix having one or more non-zero entries, whilst the n-th column in the the adjacency matrix will have only zero entries
  2. A "sink" vertex (n) will be characterised by the n-th column in the adjacency matrix having one or more non-zero entries, whilst the n-th row in the the adjacency matrix will have only zero entries
  3. If you cannot identify a "sink" and a "source" on the bove criteria, you have a graph not a network

When you have an entry on the diagonal of the adjacency matrix, this represents a "loop" and hence the associated vertex cannot be either a sink or a source, since it fails both (1) and (2) above. The mere existence of "loops" does not determine whenther or not a graph may be a network - it just identifies which verexes cannot be sinks/sources

When you are unsure whether an adjacency matrix 'A'  may represent a network or a graph it may be safer to visualise it with the DrawGraph(A) command rather than DrawNetwork(A). The former will (always?) work but the latter is pretty much guaranteed to fail unless the adjacency matrix definitely represents a network (see above criteria)

It would also be very useful if you could use the big green up-arrow in the MaplePrimes toolbar to upload Maple worksheets. A an editable/executable worksheet is just so much easier to deal with than a non-editable/non-executable, "picture" of a worksheet

Still not clear to me...

tomleslie 13876
March 11 2019
0 0

@student_md 

You can define the matrices and use these to generate the equations, or you can (somehow) enter the equations "explicitly" and from these, extract the required matrices. The attached does the following

  1. Define the matrices
  2. Generate the equations
  3. Extract the matrices from the equations and compare with (1) above
  4. Restart
  5. Define the equations more-or-less explicitly
  6. Extract the required matrices from the definition in (5) above

Something in the attached ought to be what you want

> 

  r:=3:
#
# OP uses Gamma which is protected. Override
# protection for now
#
  local Gamma:
  myDel:=(x::integer, y::integer)-> `if`(x=y,1,0):
  M:= Matrix( r,
              r,
              (i, j)-> myDel(i,j)+2*Gamma*sin(i*Pi*v__0*tau)*sin(j*Pi*v__0*tau)
            );
  C:= Matrix( r,
              r,
              (i ,j)-> 4*Gamma*Pi*v__0*j*sin(i*Pi*v__0*tau)*cos(j*Pi*v__0*tau)
            );
  K:= Matrix( r,
              r,
              (i, j)-> myDel(i,j)*i^4*Pi^2-2*Gamma(j*Pi*v__0)^2*sin(i*Pi*v__0*tau)*sin(j*Pi*v__0*tau)
            );
#
# What is tau__f?! Is 'f' meant to be defined in a
# piecewise manner?
#
  f:= Vector( r,
              i-> piecewise
                  ( tau>=0 and tau<tau__f,
                    Gamma*Pi^3*sin(i*Pi*v__0*tau),
                    0
                  )
            );
#
# What is 'p'? Use 'r for now
#
# 'D' is protected as the differetnial operator
# (which I'm guessing we might need later) so use
# 'DD' instead
#
  p:=r:
  DD:= Matrix( r,
               r,
               (i, j)-> sin(i*Pi*eta[j])
             );
  UU:=Vector( r, i->U[i](tau));
  Mphi:=Vector(r, i->phi[i](tau));
  Mphi1:= diff( Mphi, tau);
  Mphi2:= diff( Mphi1, tau);

Matrix(3, 3, {(1, 1) = 2*GAMMA*sin(Pi*`#msub(mi("v"),mi("0"))`*tau)^2+1, (1, 2) = 2*GAMMA*sin(Pi*`#msub(mi("v"),mi("0"))`*tau)*sin(2*Pi*`#msub(mi("v"),mi("0"))`*tau), (1, 3) = 2*GAMMA*sin(Pi*`#msub(mi("v"),mi("0"))`*tau)*sin(3*Pi*`#msub(mi("v"),mi("0"))`*tau), (2, 1) = 2*GAMMA*sin(Pi*`#msub(mi("v"),mi("0"))`*tau)*sin(2*Pi*`#msub(mi("v"),mi("0"))`*tau), (2, 2) = 2*GAMMA*sin(2*Pi*`#msub(mi("v"),mi("0"))`*tau)^2+1, (2, 3) = 2*GAMMA*sin(2*Pi*`#msub(mi("v"),mi("0"))`*tau)*sin(3*Pi*`#msub(mi("v"),mi("0"))`*tau), (3, 1) = 2*GAMMA*sin(Pi*`#msub(mi("v"),mi("0"))`*tau)*sin(3*Pi*`#msub(mi("v"),mi("0"))`*tau), (3, 2) = 2*GAMMA*sin(2*Pi*`#msub(mi("v"),mi("0"))`*tau)*sin(3*Pi*`#msub(mi("v"),mi("0"))`*tau), (3, 3) = 2*GAMMA*sin(3*Pi*`#msub(mi("v"),mi("0"))`*tau)^2+1})

  

Matrix(3, 3, {(1, 1) = 4*GAMMA*Pi*`#msub(mi("v"),mi("0"))`*sin(Pi*`#msub(mi("v"),mi("0"))`*tau)*cos(Pi*`#msub(mi("v"),mi("0"))`*tau), (1, 2) = 8*GAMMA*Pi*`#msub(mi("v"),mi("0"))`*sin(Pi*`#msub(mi("v"),mi("0"))`*tau)*cos(2*Pi*`#msub(mi("v"),mi("0"))`*tau), (1, 3) = 12*GAMMA*Pi*`#msub(mi("v"),mi("0"))`*sin(Pi*`#msub(mi("v"),mi("0"))`*tau)*cos(3*Pi*`#msub(mi("v"),mi("0"))`*tau), (2, 1) = 4*GAMMA*Pi*`#msub(mi("v"),mi("0"))`*sin(2*Pi*`#msub(mi("v"),mi("0"))`*tau)*cos(Pi*`#msub(mi("v"),mi("0"))`*tau), (2, 2) = 8*GAMMA*Pi*`#msub(mi("v"),mi("0"))`*sin(2*Pi*`#msub(mi("v"),mi("0"))`*tau)*cos(2*Pi*`#msub(mi("v"),mi("0"))`*tau), (2, 3) = 12*GAMMA*Pi*`#msub(mi("v"),mi("0"))`*sin(2*Pi*`#msub(mi("v"),mi("0"))`*tau)*cos(3*Pi*`#msub(mi("v"),mi("0"))`*tau), (3, 1) = 4*GAMMA*Pi*`#msub(mi("v"),mi("0"))`*sin(3*Pi*`#msub(mi("v"),mi("0"))`*tau)*cos(Pi*`#msub(mi("v"),mi("0"))`*tau), (3, 2) = 8*GAMMA*Pi*`#msub(mi("v"),mi("0"))`*sin(3*Pi*`#msub(mi("v"),mi("0"))`*tau)*cos(2*Pi*`#msub(mi("v"),mi("0"))`*tau), (3, 3) = 12*GAMMA*Pi*`#msub(mi("v"),mi("0"))`*sin(3*Pi*`#msub(mi("v"),mi("0"))`*tau)*cos(3*Pi*`#msub(mi("v"),mi("0"))`*tau)})

  

Matrix(3, 3, {(1, 1) = Pi^2-2*GAMMA(Pi*`#msub(mi("v"),mi("0"))`)^2*sin(Pi*`#msub(mi("v"),mi("0"))`*tau)^2, (1, 2) = -2*GAMMA(2*Pi*`#msub(mi("v"),mi("0"))`)^2*sin(Pi*`#msub(mi("v"),mi("0"))`*tau)*sin(2*Pi*`#msub(mi("v"),mi("0"))`*tau), (1, 3) = -2*GAMMA(3*Pi*`#msub(mi("v"),mi("0"))`)^2*sin(Pi*`#msub(mi("v"),mi("0"))`*tau)*sin(3*Pi*`#msub(mi("v"),mi("0"))`*tau), (2, 1) = -2*GAMMA(Pi*`#msub(mi("v"),mi("0"))`)^2*sin(2*Pi*`#msub(mi("v"),mi("0"))`*tau)*sin(Pi*`#msub(mi("v"),mi("0"))`*tau), (2, 2) = 16*Pi^2-2*GAMMA(2*Pi*`#msub(mi("v"),mi("0"))`)^2*sin(2*Pi*`#msub(mi("v"),mi("0"))`*tau)^2, (2, 3) = -2*GAMMA(3*Pi*`#msub(mi("v"),mi("0"))`)^2*sin(2*Pi*`#msub(mi("v"),mi("0"))`*tau)*sin(3*Pi*`#msub(mi("v"),mi("0"))`*tau), (3, 1) = -2*GAMMA(Pi*`#msub(mi("v"),mi("0"))`)^2*sin(3*Pi*`#msub(mi("v"),mi("0"))`*tau)*sin(Pi*`#msub(mi("v"),mi("0"))`*tau), (3, 2) = -2*GAMMA(2*Pi*`#msub(mi("v"),mi("0"))`)^2*sin(3*Pi*`#msub(mi("v"),mi("0"))`*tau)*sin(2*Pi*`#msub(mi("v"),mi("0"))`*tau), (3, 3) = 81*Pi^2-2*GAMMA(3*Pi*`#msub(mi("v"),mi("0"))`)^2*sin(3*Pi*`#msub(mi("v"),mi("0"))`*tau)^2})

  

Vector(3, {(1) = piecewise(0 <= tau and tau < `#msub(mi("&tau;",fontstyle = "normal"),mi("f"))`, GAMMA*Pi^3*sin(Pi*`#msub(mi("v"),mi("0"))`*tau), 0), (2) = piecewise(0 <= tau and tau < `#msub(mi("&tau;",fontstyle = "normal"),mi("f"))`, GAMMA*Pi^3*sin(2*Pi*`#msub(mi("v"),mi("0"))`*tau), 0), (3) = piecewise(0 <= tau and tau < `#msub(mi("&tau;",fontstyle = "normal"),mi("f"))`, GAMMA*Pi^3*sin(3*Pi*`#msub(mi("v"),mi("0"))`*tau), 0)})

  

Matrix(3, 3, {(1, 1) = sin(Pi*eta[1]), (1, 2) = sin(Pi*eta[2]), (1, 3) = sin(Pi*eta[3]), (2, 1) = sin(2*Pi*eta[1]), (2, 2) = sin(2*Pi*eta[2]), (2, 3) = sin(2*Pi*eta[3]), (3, 1) = sin(3*Pi*eta[1]), (3, 2) = sin(3*Pi*eta[2]), (3, 3) = sin(3*Pi*eta[3])})

  

Vector(3, {(1) = U[1](tau), (2) = U[2](tau), (3) = U[3](tau)})

  

Vector(3, {(1) = phi[1](tau), (2) = phi[2](tau), (3) = phi[3](tau)})

  

Vector(3, {(1) = diff(phi[1](tau), tau), (2) = diff(phi[2](tau), tau), (3) = diff(phi[3](tau), tau)})

  

Vector[column](%id = 18446744074584482870)

 (1)
> 

  with(LinearAlgebra):
###########################################
# Using the above matrices, produce the LHS
# of the OP's equation as a vector
#
# Objective(?!) is to extract the above matrices
# purely from this vector of equations. Since
# the required matrices are already "known",
# this seems a little pointless!!!
#
  Eqs:=M.Mphi2+C.Mphi1+K.Mphi=f+DD.UU:
#
# Define the variables
#
  vars:=[seq(phi[n](tau), n=1..r)]:
#
# From the equations Eqs extract the Matrix
# of coefficients for 'phi'
#
# Check that this is the same as 'K' defined above
#
  GenerateMatrix
  ( convert
    ( lhs(Eqs),
      list
    ),
    vars
  )[1];
  Equal(%,K);
#
# From the equations Eqs extract the Matrix
# of coefficients for 'diff(phi(tau),tau)'
#
# Check that this is the same as 'C' defined above
#
  GenerateMatrix
  ( convert
    ( lhs(Eqs),
      list
    ),
    diff(vars,tau)
  )[1];
  Equal(%,C);
#
# From the equations Eqs extract the Matrix
# of coefficients for 'diff(phi(tau),tau)'
#
# Check that this is the same as 'M' defined above
#
  GenerateMatrix
  ( convert
    ( lhs(Eqs),
      list
    ),
    diff(vars,tau,tau)
  )[1];
  Equal(%,M);
#
# Get  the terms on the rhs of the quations
#
  Vector[column]
  ( [ op~( 1,
           convert(rhs~(Eqs),list)
         )
    ]
  );
  Equal(%,f);

  GenerateMatrix
  ( convert(rhs(Eqs), list),
    [ seq(U[i](tau), i=1..r)]
  )[1];
  Equal(%,DD);

Matrix(3, 3, {(1, 1) = Pi^2-2*GAMMA(Pi*`#msub(mi("v"),mi("0"))`)^2*sin(Pi*`#msub(mi("v"),mi("0"))`*tau)^2, (1, 2) = -2*GAMMA(2*Pi*`#msub(mi("v"),mi("0"))`)^2*sin(Pi*`#msub(mi("v"),mi("0"))`*tau)*sin(2*Pi*`#msub(mi("v"),mi("0"))`*tau), (1, 3) = -2*GAMMA(3*Pi*`#msub(mi("v"),mi("0"))`)^2*sin(Pi*`#msub(mi("v"),mi("0"))`*tau)*sin(3*Pi*`#msub(mi("v"),mi("0"))`*tau), (2, 1) = -2*GAMMA(Pi*`#msub(mi("v"),mi("0"))`)^2*sin(2*Pi*`#msub(mi("v"),mi("0"))`*tau)*sin(Pi*`#msub(mi("v"),mi("0"))`*tau), (2, 2) = 16*Pi^2-2*GAMMA(2*Pi*`#msub(mi("v"),mi("0"))`)^2*sin(2*Pi*`#msub(mi("v"),mi("0"))`*tau)^2, (2, 3) = -2*GAMMA(3*Pi*`#msub(mi("v"),mi("0"))`)^2*sin(2*Pi*`#msub(mi("v"),mi("0"))`*tau)*sin(3*Pi*`#msub(mi("v"),mi("0"))`*tau), (3, 1) = -2*GAMMA(Pi*`#msub(mi("v"),mi("0"))`)^2*sin(3*Pi*`#msub(mi("v"),mi("0"))`*tau)*sin(Pi*`#msub(mi("v"),mi("0"))`*tau), (3, 2) = -2*GAMMA(2*Pi*`#msub(mi("v"),mi("0"))`)^2*sin(3*Pi*`#msub(mi("v"),mi("0"))`*tau)*sin(2*Pi*`#msub(mi("v"),mi("0"))`*tau), (3, 3) = 81*Pi^2-2*GAMMA(3*Pi*`#msub(mi("v"),mi("0"))`)^2*sin(3*Pi*`#msub(mi("v"),mi("0"))`*tau)^2})

  

true

  

Matrix(3, 3, {(1, 1) = 4*GAMMA*Pi*`#msub(mi("v"),mi("0"))`*sin(Pi*`#msub(mi("v"),mi("0"))`*tau)*cos(Pi*`#msub(mi("v"),mi("0"))`*tau), (1, 2) = 8*GAMMA*Pi*`#msub(mi("v"),mi("0"))`*sin(Pi*`#msub(mi("v"),mi("0"))`*tau)*cos(2*Pi*`#msub(mi("v"),mi("0"))`*tau), (1, 3) = 12*GAMMA*Pi*`#msub(mi("v"),mi("0"))`*sin(Pi*`#msub(mi("v"),mi("0"))`*tau)*cos(3*Pi*`#msub(mi("v"),mi("0"))`*tau), (2, 1) = 4*GAMMA*Pi*`#msub(mi("v"),mi("0"))`*sin(2*Pi*`#msub(mi("v"),mi("0"))`*tau)*cos(Pi*`#msub(mi("v"),mi("0"))`*tau), (2, 2) = 8*GAMMA*Pi*`#msub(mi("v"),mi("0"))`*sin(2*Pi*`#msub(mi("v"),mi("0"))`*tau)*cos(2*Pi*`#msub(mi("v"),mi("0"))`*tau), (2, 3) = 12*GAMMA*Pi*`#msub(mi("v"),mi("0"))`*sin(2*Pi*`#msub(mi("v"),mi("0"))`*tau)*cos(3*Pi*`#msub(mi("v"),mi("0"))`*tau), (3, 1) = 4*GAMMA*Pi*`#msub(mi("v"),mi("0"))`*sin(3*Pi*`#msub(mi("v"),mi("0"))`*tau)*cos(Pi*`#msub(mi("v"),mi("0"))`*tau), (3, 2) = 8*GAMMA*Pi*`#msub(mi("v"),mi("0"))`*sin(3*Pi*`#msub(mi("v"),mi("0"))`*tau)*cos(2*Pi*`#msub(mi("v"),mi("0"))`*tau), (3, 3) = 12*GAMMA*Pi*`#msub(mi("v"),mi("0"))`*sin(3*Pi*`#msub(mi("v"),mi("0"))`*tau)*cos(3*Pi*`#msub(mi("v"),mi("0"))`*tau)})

  

true

  

Matrix(3, 3, {(1, 1) = 2*GAMMA*sin(Pi*`#msub(mi("v"),mi("0"))`*tau)^2+1, (1, 2) = 2*GAMMA*sin(Pi*`#msub(mi("v"),mi("0"))`*tau)*sin(2*Pi*`#msub(mi("v"),mi("0"))`*tau), (1, 3) = 2*GAMMA*sin(Pi*`#msub(mi("v"),mi("0"))`*tau)*sin(3*Pi*`#msub(mi("v"),mi("0"))`*tau), (2, 1) = 2*GAMMA*sin(Pi*`#msub(mi("v"),mi("0"))`*tau)*sin(2*Pi*`#msub(mi("v"),mi("0"))`*tau), (2, 2) = 2*GAMMA*sin(2*Pi*`#msub(mi("v"),mi("0"))`*tau)^2+1, (2, 3) = 2*GAMMA*sin(2*Pi*`#msub(mi("v"),mi("0"))`*tau)*sin(3*Pi*`#msub(mi("v"),mi("0"))`*tau), (3, 1) = 2*GAMMA*sin(Pi*`#msub(mi("v"),mi("0"))`*tau)*sin(3*Pi*`#msub(mi("v"),mi("0"))`*tau), (3, 2) = 2*GAMMA*sin(2*Pi*`#msub(mi("v"),mi("0"))`*tau)*sin(3*Pi*`#msub(mi("v"),mi("0"))`*tau), (3, 3) = 2*GAMMA*sin(3*Pi*`#msub(mi("v"),mi("0"))`*tau)^2+1})

  

true

  

Vector(3, {(1) = piecewise(0 <= tau and tau < `#msub(mi("&tau;",fontstyle = "normal"),mi("f"))`, GAMMA*Pi^3*sin(Pi*`#msub(mi("v"),mi("0"))`*tau), 0), (2) = piecewise(0 <= tau and tau < `#msub(mi("&tau;",fontstyle = "normal"),mi("f"))`, GAMMA*Pi^3*sin(2*Pi*`#msub(mi("v"),mi("0"))`*tau), 0), (3) = piecewise(0 <= tau and tau < `#msub(mi("&tau;",fontstyle = "normal"),mi("f"))`, GAMMA*Pi^3*sin(3*Pi*`#msub(mi("v"),mi("0"))`*tau), 0)})

  

true

  

Matrix(%id = 18446744074583731246)

  

true

 (2)
> 

#######################################################
# Alternatively one can write the vector of equations
# more or less explicitly and then extract the
# required matrices
#
  restart;
  with(LinearAlgebra):
  r:=3:
  local Gamma:
#
# Write the vector of equations
#
  Eqs:= Vector[column]
        ( r,
          i-> diff(phi[i](tau), tau,tau)
              +
              2*Gamma*add( sin(i*Pi*v__0*tau)*sin(j*Pi*v__0*tau)*diff(phi[j](tau), tau,tau),
                           j=1..r
                         )
              +
              4*Gamma*Pi*v__0*add( j*sin(i*Pi*v__0*tau)*cos(j*Pi*v__0*tau)*diff(phi[j](tau),tau),

                                   j=1..r
                                 )
              +
              i^4*Pi^2*phi[i](tau)
              -
              2*Gamma*add( (j*Pi*v__0)^3*sin(i*Pi*v__0*tau)*sin(j*Pi*v__0*tau)*phi[j](tau),
                           j=1..r
                         )
              =
              piecewise( tau>0 and tau<tau__f,
                         (Gamma*Pi^3/2)*sin(i*Pi*v__0*tau),
                         0
                       )
              +
              add( U[j](tau)*sin(i*Pi*eta[j]), j=1..r)
        ):
#
# Now extract the various matrices
#
  M:= GenerateMatrix
      ( lhs~(convert(Eqs, list)),
        [seq(phi[n](tau), n=1..r)]
      )[1];
  C:= GenerateMatrix
      ( lhs~(convert(Eqs, list)),
        diff([seq(phi[n](tau), n=1..r)],tau)
      )[1];
  K:= GenerateMatrix
      ( lhs~(convert(Eqs, list)),
        diff([seq(phi[n](tau), n=1..r)],tau,tau)
      )[1];
  f:= Vector[column]
      ( [ op~( 1,
               convert(rhs~(Eqs),list)
             )
        ]
      );
  DD:= GenerateMatrix
       ( rhs~(convert(Eqs, list)),
         [seq(U[i](tau), i=1..r)]
       )[1];
  

Matrix(3, 3, {(1, 1) = -2*GAMMA*sin(Pi*`#msub(mi("v"),mi("0"))`*tau)^2*Pi^3*`#msub(mi("v"),mi("0"))`^3+Pi^2, (1, 2) = -16*GAMMA*sin(Pi*`#msub(mi("v"),mi("0"))`*tau)*sin(2*Pi*`#msub(mi("v"),mi("0"))`*tau)*Pi^3*`#msub(mi("v"),mi("0"))`^3, (1, 3) = -54*GAMMA*Pi^3*`#msub(mi("v"),mi("0"))`^3*sin(Pi*`#msub(mi("v"),mi("0"))`*tau)*sin(3*Pi*`#msub(mi("v"),mi("0"))`*tau), (2, 1) = -2*GAMMA*sin(Pi*`#msub(mi("v"),mi("0"))`*tau)*sin(2*Pi*`#msub(mi("v"),mi("0"))`*tau)*Pi^3*`#msub(mi("v"),mi("0"))`^3, (2, 2) = -16*GAMMA*sin(2*Pi*`#msub(mi("v"),mi("0"))`*tau)^2*Pi^3*`#msub(mi("v"),mi("0"))`^3+16*Pi^2, (2, 3) = -54*GAMMA*Pi^3*`#msub(mi("v"),mi("0"))`^3*sin(2*Pi*`#msub(mi("v"),mi("0"))`*tau)*sin(3*Pi*`#msub(mi("v"),mi("0"))`*tau), (3, 1) = -2*GAMMA*Pi^3*`#msub(mi("v"),mi("0"))`^3*sin(Pi*`#msub(mi("v"),mi("0"))`*tau)*sin(3*Pi*`#msub(mi("v"),mi("0"))`*tau), (3, 2) = -16*GAMMA*Pi^3*`#msub(mi("v"),mi("0"))`^3*sin(2*Pi*`#msub(mi("v"),mi("0"))`*tau)*sin(3*Pi*`#msub(mi("v"),mi("0"))`*tau), (3, 3) = -54*GAMMA*Pi^3*`#msub(mi("v"),mi("0"))`^3*sin(3*Pi*`#msub(mi("v"),mi("0"))`*tau)^2+81*Pi^2})

  

Matrix(3, 3, {(1, 1) = 4*GAMMA*Pi*`#msub(mi("v"),mi("0"))`*sin(Pi*`#msub(mi("v"),mi("0"))`*tau)*cos(Pi*`#msub(mi("v"),mi("0"))`*tau), (1, 2) = 8*GAMMA*Pi*`#msub(mi("v"),mi("0"))`*sin(Pi*`#msub(mi("v"),mi("0"))`*tau)*cos(2*Pi*`#msub(mi("v"),mi("0"))`*tau), (1, 3) = 12*GAMMA*Pi*`#msub(mi("v"),mi("0"))`*sin(Pi*`#msub(mi("v"),mi("0"))`*tau)*cos(3*Pi*`#msub(mi("v"),mi("0"))`*tau), (2, 1) = 4*GAMMA*Pi*`#msub(mi("v"),mi("0"))`*sin(2*Pi*`#msub(mi("v"),mi("0"))`*tau)*cos(Pi*`#msub(mi("v"),mi("0"))`*tau), (2, 2) = 8*GAMMA*Pi*`#msub(mi("v"),mi("0"))`*sin(2*Pi*`#msub(mi("v"),mi("0"))`*tau)*cos(2*Pi*`#msub(mi("v"),mi("0"))`*tau), (2, 3) = 12*GAMMA*Pi*`#msub(mi("v"),mi("0"))`*sin(2*Pi*`#msub(mi("v"),mi("0"))`*tau)*cos(3*Pi*`#msub(mi("v"),mi("0"))`*tau), (3, 1) = 4*GAMMA*Pi*`#msub(mi("v"),mi("0"))`*sin(3*Pi*`#msub(mi("v"),mi("0"))`*tau)*cos(Pi*`#msub(mi("v"),mi("0"))`*tau), (3, 2) = 8*GAMMA*Pi*`#msub(mi("v"),mi("0"))`*sin(3*Pi*`#msub(mi("v"),mi("0"))`*tau)*cos(2*Pi*`#msub(mi("v"),mi("0"))`*tau), (3, 3) = 12*GAMMA*Pi*`#msub(mi("v"),mi("0"))`*sin(3*Pi*`#msub(mi("v"),mi("0"))`*tau)*cos(3*Pi*`#msub(mi("v"),mi("0"))`*tau)})

  

Matrix(3, 3, {(1, 1) = 2*GAMMA*sin(Pi*`#msub(mi("v"),mi("0"))`*tau)^2+1, (1, 2) = 2*GAMMA*sin(Pi*`#msub(mi("v"),mi("0"))`*tau)*sin(2*Pi*`#msub(mi("v"),mi("0"))`*tau), (1, 3) = 2*GAMMA*sin(Pi*`#msub(mi("v"),mi("0"))`*tau)*sin(3*Pi*`#msub(mi("v"),mi("0"))`*tau), (2, 1) = 2*GAMMA*sin(Pi*`#msub(mi("v"),mi("0"))`*tau)*sin(2*Pi*`#msub(mi("v"),mi("0"))`*tau), (2, 2) = 2*GAMMA*sin(2*Pi*`#msub(mi("v"),mi("0"))`*tau)^2+1, (2, 3) = 2*GAMMA*sin(2*Pi*`#msub(mi("v"),mi("0"))`*tau)*sin(3*Pi*`#msub(mi("v"),mi("0"))`*tau), (3, 1) = 2*GAMMA*sin(Pi*`#msub(mi("v"),mi("0"))`*tau)*sin(3*Pi*`#msub(mi("v"),mi("0"))`*tau), (3, 2) = 2*GAMMA*sin(2*Pi*`#msub(mi("v"),mi("0"))`*tau)*sin(3*Pi*`#msub(mi("v"),mi("0"))`*tau), (3, 3) = 2*GAMMA*sin(3*Pi*`#msub(mi("v"),mi("0"))`*tau)^2+1})

  

Vector(3, {(1) = piecewise(0 < tau and tau < `#msub(mi("&tau;",fontstyle = "normal"),mi("f"))`, (1/2)*GAMMA*Pi^3*sin(Pi*`#msub(mi("v"),mi("0"))`*tau), 0), (2) = piecewise(0 < tau and tau < `#msub(mi("&tau;",fontstyle = "normal"),mi("f"))`, (1/2)*GAMMA*Pi^3*sin(2*Pi*`#msub(mi("v"),mi("0"))`*tau), 0), (3) = piecewise(0 < tau and tau < `#msub(mi("&tau;",fontstyle = "normal"),mi("f"))`, (1/2)*GAMMA*Pi^3*sin(3*Pi*`#msub(mi("v"),mi("0"))`*tau), 0)})

  

Matrix(%id = 18446744074408680374)

 (3)
> 

 

Download setup2.mw

Observation...

tomleslie 13876
March 10 2019
0 0

@mmcdara 

As supplied, the matrix A defines a directed graph - for example A[2,4] is not equal to A[4,2], hence the adjacency matrix is not symmetric, hence the graph is directed.

Observation...

tomleslie 13876
March 09 2019
0 0

@acer 
The Support() command does not only "serve(s) for the uniform case". It works for all distributions and returns the range of values for which the probability density function is non-zero. Suggest you read the help.

Even for continuous (non-uniform) distributions, it can be infomative. EG fo the "normal" distribution

with(Statistics):
Support(Normal(3, 1), 'output = range');

it will correctly return -inifinity..infinity.

My only reason for using the Support() function in my response was that it was a convenient way to handle the issue of dice which had <>6 faces. This is not "abstruse" - just convenient

 

 

I have solved...

tomleslie 13876
March 07 2019
0 0

@Gabriel samaila 

the problem which you have posted

It is surely obvious that I cannot solve some problem which you have not posted

Not how things work...

tomleslie 13876
March 07 2019
0 0

@Gabriel samaila 

You are telling me that you do not even know what problem you are trying to solve!!!

If you don't know - how do you expect me to know???!!!

Maple's "automatic" simplification...

tomleslie 13876
March 05 2019
0 0

@erik10 

is sometimes a benefit and sometimes a real pain. It is easy to forget that it is pretty much always happening in the background. It is certainky dangerous to rely on it in any code. Apart from anything else, I don't think that exactly what is "simplified" is listed/described anywhere in the documentation. Just for your amusement, contemplate the expression

DIce1()/Dice1();

If this is "simplified" before any of the calls to Dice1() are "evaluated", then the answer will always be 1. If it is not "simplified", then you will get all combinations of rationals from 1/6 to 6. In fact the latter happens.

So Dice1()+Dice1() "simplifies" to 2*Dice1(), but Dice1()/Dice1() does not "simplify" to 1.

Given that Dice1()+Dice1() "simplifies" would you like to place a bet on the outcome of Dice1()-Dice1() without actually trying it???

First 71 72 73 74 75 76 77 Last Page 73 of 207