0

I am having some problems in plotting this. Everything is ok until the plot statement where nothing plots. Can someone please help me so that it can plot something. The following is my code:

j = 10; 
s = 0; r = 0;

B[n_] = Integrate[2*Sin[n*Pi*x]*(x), {x, 0, 1}];
u[x_, psi_] = Sum[B[n]*Sin[n*Pi*x]*Exp[-(n*Pi)^2*psi], {n, 1, j}];

K[x_, psi_] = 
  Sum[Sin[n*Pi*x]*
    Sin[n*Pi*
      psi]*(2*Exp[-(n*Pi)^2*
         Abs[s + r]] - (Exp[-(n*Pi)^2*Abs[s - r]] - 
         Exp[-(n*Pi)^2*(s + r)])/(n*Pi)^2 ), {n, 1, j}];

w = RandomReal[NormalDistribution[0, 1], 101];
d = Round[100*x + 1];
S = Total[Total[u[x, psi]/Length[u[x, psi]]] + w[d]]

T[x_, psi_] = Integrate[K[x - y, psi]*(y)*S, {y, -10, 10}]

Plot3D[T[x, psi], {x, 0, 1}, {psi, 0.01, 1}, 
 AxesLabel -> {"x", "t", "Temperature"}, Boxed -> False, 
 Mesh -> False]

Basically, I have some data from "u" and I want to make it noisy (from "w") for each "x" value and then perform the convolution in "T" and plot.

I will really appreciate anyone's kind help.

Thanks very much!

user1482746
  • 5
  • 1
  • 4
  • In the definition of `S` you write `w[d]` but `w` is a list of numbers, not a function. You could try something like `Interpolation[w][d]` instead. – Heike Jul 11 '12 at 16:47

3 Answers3

1

It looks very much as if you are using = where you should be using :=. The former makes an immediate assignment (called Set) the other a delayed assignment (SetDelayed). The difference is fundamental in Mathematica, you should read the documentation until you understand this difference thoroughly.

High Performance Mark
  • 77,191
  • 7
  • 105
  • 161
  • can you please tell me where I should use := for it to work. Sorry about the inconvenience. – user1482746 Jul 11 '12 at 14:39
  • I have been getting the following error: Part::pspec saying it is not an integer in "S". The integrated solution of "T" looks fine as well, but I have no clue why the plot doesn't work because it makes sense as plugging in numbers and plotting. – user1482746 Jul 11 '12 at 14:45
  • Have you seen `Convolve` and `ListConvolve`? – image_doctor Jul 11 '12 at 15:11
  • I did try to use that first, but it didn't work for me so I did it this way. If anyone can please try to correct me where I may be going wrong. I'll greatly appreciate it. – user1482746 Jul 11 '12 at 15:36
  • A link to the documentation which pertains to this topic would improve this answer. – Justin Oroz Mar 24 '17 at 06:21
1

I'm not sure that I understand the problem you're trying to solve. However, modifying your code as shown below allows it to run - I rephrased several expression to be functions (a good rule of thumb is to use := if the left hand side involves a pattern, like B[n_]) and I removed some code that was apparently trying to treat scalars as vectors.

j = 10; s = 0; r = 0;

ClearAll[B];
B[n_] := B[n] = Integrate[2*Sin[n*Pi*a]*(a), {a, 0, 1}];

ClearAll[u];
u[x_, psi_] := Sum[B[n]*Sin[n*Pi*x]*Exp[-(n*Pi)^2*psi], {n, 1, j}];

K[x_, psi_] := 
  Sum[Sin[n*Pi*x]*
    Sin[n*Pi*
      psi]*(2*Exp[-(n*Pi)^2*Abs[s + r]] - (Exp[-(n*Pi)^2*Abs[s - r]] -
          Exp[-(n*Pi)^2*(s + r)])/(n*Pi)^2), {n, 1, j}];

S[x_, psi_] := u[x, psi] + RandomReal[NormalDistribution[]]

T[x_, psi_] := Integrate[K[x - y, psi]*(y)*S[x, psi], {y, -10, 10}]

Plot3D[T[x, psi], {x, 0, 1}, {psi, 0.01, 1}, 
 AxesLabel -> {"x", "t", "Temperature"}, Boxed -> False, 
 Mesh -> False]

After running for some time (~ 1 hour) it produces the plot below

Plot of T

There is probably a much more efficient way to produce this plot using a more idiomatic approach. If you could provide more detailed information about what you're trying to do with the code you posted, then maybe I or others could give you a more useful answer.

DGrady
  • 1,065
  • 1
  • 13
  • 23
  • I tried running the code on my machine and it did not give me a plot :( just empty axes...can you please tell me as to what may be wrong? I used your modified code – user1482746 Jul 12 '12 at 18:42
  • What version of MM are you using? I'm on 8.0.1.0 running Mac OS X 10.7.4. Do you see a plot if you run Plot3D[Sin[x + y^2], {x, -3, 3}, {y, -2, 2}] ? – DGrady Jul 12 '12 at 19:18
  • yup i got that plot. I'm running MM 7.0.0 and I'm running it in Windows 7 (bootcamp of MAC) – user1482746 Jul 12 '12 at 19:26
  • `RandomVariate` was introduced in MM 8; maybe that's causing the problem? I updated the answer to use RandomReal instead. – DGrady Jul 12 '12 at 19:38
  • thanks soo much...it works! I was also wondering if you have some time to space and if you could look a the other question I posted, sorry to bother you even though you helped so much already. I'll deeply appreciate it: http://mathematica.stackexchange.com/questions/8251/mathematica-plotting-of-complex-and-convolution-integrals – user1482746 Jul 12 '12 at 19:56
0

Here is a template solution based on the outline of your question:

data = RandomInteger[{0, 1}, 100]; (* data creation function *) 
noise = RandomVariate[NormalDistribution[0, 1], Length@data]; (* noise vector *)
noisyData = data + noise; (* sum noise and data *)
ListConvolve[data, noisyData] (* apply convolution *)

{8.20928}

How does this prototype match with your goals ?

image_doctor
  • 501
  • 3
  • 6