Heat flow in 1d


I want to Plot Theta vs index " i ". where θ[i-1]+b*θ[i]+θ[i+1] =0 is the equation that i obtained after
discretization . b is parameter. i took m=100 and delta g= 0.01

using Plots

begin
	Δg=0.01 ;m=100
	b=-2-(m^2)*(Δg)^2
end

θ=zeros(Int8,100)

begin
	θ[1]=0
	θ[100]=1
end

p=plot()

for i in 2:99
	θ[i-1]+b*θ[i]+θ[i+1] =0
	print(θ[i-1])
	plot(θ[:])
end

current()

Hi @raman_kumar

There are two corrections i could find (i might be wrong though) →

  1. In screenshot it is mentioned that θᵢ= 1 for i = 1 and θᵢ = 0 for i = n
    So we should get below code right ? →

    begin
    	θ[1] = 1
    	θ[100] = 0
    end
    
  2. θ[i-1] + b*θ[i] + θ[i+1] = 0 is not a valid julia syntax, if you need to update θ array by indexing, probably we need to do as follows ? →

    for i in 2:99
    	θ[i] = (-1/b) * (θ[i-1] + θ[i+1])
    	print(θ[i-1], " ")
    	plot(θ[:])
    end
    

These are just my observations !!

Values of θ should be only Int8 ? because if we change the array to hold Float64, we get the values from θ[2] = 0.3333333333333333 to θ[99] = 1.7462926957343605e-47 with the same graph as you posted (decreasing trend)

image
i am not able to solve round error.

If you change the θ array to hold Float64 values, you won’t get this error…
From →

θ=zeros(Int8,100)

to →

θ=zeros(Float64,100)  # or θ=zeros(100), by default we get Float64 values

Currently it is trying to convert final calculated values (which are floats…) into Int8 (bcz your array elements are of type Int8) values which is not possible, so by changing to Float64 - θ values you would get correct result without type conversions or rounding.

1 Like


yes i have solved that issue of Float64. but still there is change in graph as i run the last code again and again as values of theta are increasing every time i want the plot at the end so that it looks exponentially decreasing.

could you please paste the code as well, so that it is easy to copy paste and try,
u can do so inside start and end triple back ticks (```) like shown below →
```
insert your code here
```

using Plots

begin
Δg=0.000000001 ;m=90
b=-2-(m^2)*(Δg)^2
end

θ=zeros(100)

begin
θ[1]=1;θ[100]=0
end

p=plot();

for i in 2:99

	while i<99
		θ[i] = round((-1/b)*(θ[i-1] + θ[i+1]),digits=9);
        print(θ[i]);
            i+=10^60
        plot(θ[:], seriestype=:scatter)
	end
end
  • when you will run last code again and again then the values of theta will increase. *
1 Like

Yes, you are right, because θ is an array which is mutable (we can modify its contents) so when you run the code for 2nd time, it would use the θ values from 1st iteration (which are not zeros anymore) and so on…
so probably u just need to update the whole array once in one loop and then you can plot the complete θ values just once, which would look like exponential decreasing scatter plot.
Code is as follows →

using Plots

begin
Δg=0.000000001 ;m=90
b=-2-(m^2)*(Δg)^2
end

θ=zeros(100)

begin
θ[1]=1;θ[100]=0
end

p=plot();

for i in 2:99
    θ[i] = round((-1/b)*(θ[i-1] + θ[i+1]),digits=9);
    print(θ[i]);
end

plot(θ[:], seriestype=:scatter)

** That Plot is fine and beautiful. Thank you Karthik **

1 Like

You are welcome :slight_smile: