Carl Love

Carl Love

26513 Reputation

25 Badges

11 years, 190 days
Himself
Wayland, Massachusetts, United States
My name was formerly Carl Devore.

MaplePrimes Activity


These are replies submitted by Carl Love

Would it be possible to give a numeric value for the parameter g?

You mentioned wanting to do this with a seq command. This is possible, even in one line. This method is even significantly faster than my previous answer. But this method can produce a tremendous amount of garbage (memory) to be collected, because it generates the powerset of T. So, WARNING: Do not use this procedure if is larger than about 25 sets.

UnionClosure:= (T::set(set))->
   {seq}(`union`(X[]), X= combinat:-powerset(T minus {{}, `union`(T[])})):


UnionClosure({{1}, {1,2}, {3}});
           {{}, {1}, {3}, {1, 2}, {1, 3}, {1, 2, 3}}

Both of my UnionClosure procedures are many, many times more efficient than the previous answers in this thread.

You mentioned wanting to do this with a seq command. This is possible, even in one line. This method is even significantly faster than my previous answer. But this method can produce a tremendous amount of garbage (memory) to be collected, because it generates the powerset of T. So, WARNING: Do not use this procedure if is larger than about 25 sets.

UnionClosure:= (T::set(set))->
   {seq}(`union`(X[]), X= combinat:-powerset(T minus {{}, `union`(T[])})):


UnionClosure({{1}, {1,2}, {3}});
           {{}, {1}, {3}, {1, 2}, {1, 3}, {1, 2, 3}}

Both of my UnionClosure procedures are many, many times more efficient than the previous answers in this thread.

Throughout this Comment, let |V| represent the number of vertices in a (finite simple) graph, let |E| represent the number of edges, and let P(C) represent the probability that a random graph (randomly selected by the "G(n, p)" model of Erdos and Renyi) with those parameters is connected. As Markiyan's title of the Post suggests, the distribution of P(C) as a function of |E| for given |V| is an all-or-nothing phenomenon. Or, in the language of the Wikipedia page that was referenced, it is a "sharp threshold" phenomenon, with the threshold being |E| = ln(|V|)*(|V|-1)/2. (This corresponds to n = |V|, p = ln(|V|) / |V| in the G(n, p) model.) In other words, for |E| > ln(|V|)*(|V|-1)/2, limit(P(C), |V|= infinity) = 1; and for |E| < ln(|V|)*(|V|-1)/2, limit(P(C), |V|= infinity) = 0. Or in the language of Random Graph Theory, in the former case "almost every" graph in G(n, p) is connected; while in the latter case "almost every" graph is disconnected.

In the plot, I would've liked to see the numbers of vertices and edges on the axes. In the code below I regenerate the plot with this information. Otherwise, the plot uses the same parameters as the original. I have also plotted the threshold function ln(|V|)*(|V|-1)/2 on the same plot. You can see that the plot models the threshold phenomenon nearly perfectly.

 

restart;

Ntrials:= 100:  
Emin:= 60:  Emax:= 360:  Eskip:= 30:
Vmin:= 30:  Vmax:= 360:  Vskip:= 10:

Pr:= Array(
     1 .. 1 + iquo(Vmax-Vmin, Vskip),
     1 .. 1 + iquo(Emax-Emin, Eskip),
     1..3,
     datatype= float[8]
):

for E from Emin by Eskip to Emax do
     for V from Vmin by Vskip to Vmax do
          ct:= 0;
          if E >= V-1 then  #threshold of possibility of connectedness
               to Ntrials do
                    if GraphTheory:-IsConnected(
                         GraphTheory:-RandomGraphs:-RandomGraph(V,E)
                    ) then
                         ct:= ct+1
                    end if
               end do
          end if;
          Pr[1 + (V-Vmin)/Vskip, 1 + (E-Emin)/Eskip, ..]:=

               < V, E, ct/Ntrials >
     end do
end do:

P1:= plots:-surfdata(Pr, transparency= .4):

P2:= plots:-spacecurve([n, ln(n)*(n-1)/2, 0], n= Vmin..Vmax,
     thickness= 4, color= red
):

plots:-display([P1,P2], view= [Vmin..Vmax, Emin..Emax, 0..1],
     labels= ['V','E','P(C)']
);

 

 

Download RandomGraphs.mw

@J4James Good question. Ordinarily, I wouldn't rely on "magic numbers" such as the position numbers 3 or 4. And there's no need to if you use eval. Let's take a look at the structure of the "listprocedure":

PDE:={diff(u(x,t),t)=w(x,t), diff(u(x,t),x)=-w(x,t)}:
IBC:= {u(x,0)=sin(2*Pi*x),u(0,t)=-sin(2*Pi*t)}:
pds:= pdsolve(PDE, IBC, numeric, time=t, range=0..1):

LP:= pds:-value(output=listprocedure);
  LP:= [x = proc()  ...  end proc, t = proc()  ...  end proc,

    u(x, t) = proc()  ...  end proc, w(x, t) = proc()  ...  end proc]

So, it's a list of assignment equations. So we apply eval:

W:= eval(w(x,t), LP);
W:= proc()  ...  end proc

Is that clear?

@J4James Good question. Ordinarily, I wouldn't rely on "magic numbers" such as the position numbers 3 or 4. And there's no need to if you use eval. Let's take a look at the structure of the "listprocedure":

PDE:={diff(u(x,t),t)=w(x,t), diff(u(x,t),x)=-w(x,t)}:
IBC:= {u(x,0)=sin(2*Pi*x),u(0,t)=-sin(2*Pi*t)}:
pds:= pdsolve(PDE, IBC, numeric, time=t, range=0..1):

LP:= pds:-value(output=listprocedure);
  LP:= [x = proc()  ...  end proc, t = proc()  ...  end proc,

    u(x, t) = proc()  ...  end proc, w(x, t) = proc()  ...  end proc]

So, it's a list of assignment equations. So we apply eval:

W:= eval(w(x,t), LP);
W:= proc()  ...  end proc

Is that clear?

I just want to state more clearly what taylor and numapprox:-laurent do, because it seems a bit strange (that they don't do much). Neither of these commands will ever return a series which is different from the one returned by the basic command series. They take the series returned by series, decide whether that series is a Taylor or Laurent series, and return either that series unchanged or an error message. So it is essentially a type checking of the output of series. This might be useful when used in a program, but it doesn't seem worthwhile in "desktop calculator mode" since it is usually obvious to the "naked eye" whether a series is Laurent, Taylor, or neither.

I just want to state more clearly what taylor and numapprox:-laurent do, because it seems a bit strange (that they don't do much). Neither of these commands will ever return a series which is different from the one returned by the basic command series. They take the series returned by series, decide whether that series is a Taylor or Laurent series, and return either that series unchanged or an error message. So it is essentially a type checking of the output of series. This might be useful when used in a program, but it doesn't seem worthwhile in "desktop calculator mode" since it is usually obvious to the "naked eye" whether a series is Laurent, Taylor, or neither.

Extending that technique, we get the following options:

ex:= 1/sqrt(3);
                            1  (1/2)
                            - 3     
                            3       

subs(3= `3`, 1/3= 1/`3`, ex);
                               1   
                             ------
                              (1/2)
                             3     

subs(3= ``(3), 1/3= 1/``(3), ex);
                               1    
                            --------
                               (1/2)
                            (3)     

which all appear with the standard square root sign instead of a fractional exponent in the GUI.

Extending that technique, we get the following options:

ex:= 1/sqrt(3);
                            1  (1/2)
                            - 3     
                            3       

subs(3= `3`, 1/3= 1/`3`, ex);
                               1   
                             ------
                              (1/2)
                             3     

subs(3= ``(3), 1/3= 1/``(3), ex);
                               1    
                            --------
                               (1/2)
                            (3)     

which all appear with the standard square root sign instead of a fractional exponent in the GUI.

@rachie90 The situation which is described in the paragraph following the one containing equation (2.2) in Tsuchiya's paper is exactly what is happening in your case: G'/G'' is going to 0/0 (where I'm abusing the derivative notation a bit, but you know what I mean). G'' is going to 0 a little faster, so this is manifested in the floating-point arithmetic as the values getting extremely large.

So, I guess that it's time to try equation (2.3). I think (although I'm far from sure) that this is equivalent to handling r and z together, so we're just dealing with gradients and Hessians of functions of two variables.

I guess that you saw that Tsuchiya does the catenoid problem later in the paper.

@rlopez  Thanks, it works. A plot shows clearly that -W is the correct derivative. Either

plot(-W(0,t), t= 0..1);

or

pds:-plot(-w(x,t), x=0, t= 0..1);

gives the plot of diff(u(x,t), x) for x=0.

With pdsolve, there seems to be a very subtle difference between what works in this context and what is "too far from a standard form".

@rlopez  Thanks, it works. A plot shows clearly that -W is the correct derivative. Either

plot(-W(0,t), t= 0..1);

or

pds:-plot(-w(x,t), x=0, t= 0..1);

gives the plot of diff(u(x,t), x) for x=0.

With pdsolve, there seems to be a very subtle difference between what works in this context and what is "too far from a standard form".

Thank you for your Answer. I was not familiar with the method of extracting the solution procedures. To answer the question that you pose rhetorically in the title of your Answer, there are two reasons:

  1. It just seems more natural to me to work with the plots, because I am well acquainted with the internal structure of plots and rtables, but not at all well acquainted with the internal structure of the modules returned by pdsolve. Certainly this is just a matter of taste.
  2. It hardly matters for this small example, but the plot method seems to use substantially less cputime (about 1/3 the time in this example). I haven't looked into why this is, but my guess is that the plots allow for more processing by evalhf in a "batch mode", i.e., fewer total calls to evalhf with each doing more work.

Robert Lopez wrote:

To get the values of the highest-ordered derivative, ... trick the system into thinking ... that the derivative is actually one of the system's dependent variables.

Could you give an example of this with pdsolve(..., numeric)? I am familiar with how to do this with dsolve(..., numeric), but when I try with pdsolve, I get an error message about the system being "too far from a standard form".

PDE:= {diff(u(x,t),t) = -diff(u(x,t),x), diff(u(x,t),x)= w(x,t)}:
IBC:= {u(x,0)=sin(2*Pi*x),u(0,t)=-sin(2*Pi*t)}:
pds:= pdsolve(PDE, IBC, numeric, time=t, range=0..1);
Error, (in pdsolve/numeric/par_hyp) input system is too far from a 'standard' form (see ?pdsolve,numeric for more detail)

Thank you for your Answer. I was not familiar with the method of extracting the solution procedures. To answer the question that you pose rhetorically in the title of your Answer, there are two reasons:

  1. It just seems more natural to me to work with the plots, because I am well acquainted with the internal structure of plots and rtables, but not at all well acquainted with the internal structure of the modules returned by pdsolve. Certainly this is just a matter of taste.
  2. It hardly matters for this small example, but the plot method seems to use substantially less cputime (about 1/3 the time in this example). I haven't looked into why this is, but my guess is that the plots allow for more processing by evalhf in a "batch mode", i.e., fewer total calls to evalhf with each doing more work.

Robert Lopez wrote:

To get the values of the highest-ordered derivative, ... trick the system into thinking ... that the derivative is actually one of the system's dependent variables.

Could you give an example of this with pdsolve(..., numeric)? I am familiar with how to do this with dsolve(..., numeric), but when I try with pdsolve, I get an error message about the system being "too far from a standard form".

PDE:= {diff(u(x,t),t) = -diff(u(x,t),x), diff(u(x,t),x)= w(x,t)}:
IBC:= {u(x,0)=sin(2*Pi*x),u(0,t)=-sin(2*Pi*t)}:
pds:= pdsolve(PDE, IBC, numeric, time=t, range=0..1);
Error, (in pdsolve/numeric/par_hyp) input system is too far from a 'standard' form (see ?pdsolve,numeric for more detail)

First 642 643 644 645 646 647 648 Last Page 644 of 689