On Dec 15, 2008, at 5:20 PM, Juliane Struve wrote:
> Dear list,
>
> I am trying to program semi-random movement within a circle, with no
> particles leaving the circle. I would like them to bounce back when
> they come to close to the wall, but I don't seem to be able to get
> this right. Would somebody be able to give me a hint ? This is my
> code so far, the particle starts at some point and moves towards the
> wall, but I don't get the "bouncing off" part right .
That is a bit vague for an audience that is numerically oriented.
>
>
> Any help would be much appreciated.
>
> Juliane
>
> days=10
> circularspace
> =
> data
> .frame
> (day
> =c(0:days),xcoord=1,ycoord=1,xvelocity=1,yvelocity=1,xdistwall=0,
> ydistwall=0, wallxvel=0, wallyvel=0,stochasticxvel=0,stochasticyvel=0)
You have initialized this object with 10 rows, but why would you
specify the distance to the walls as 0?
>
> xmax=10
> xmin=-10
> ymax=10
> ymin=-10
> mindist=8
> plot(xmin:xmax, ymin:ymax, type = "n")
> circularspace
> radius=10
> timesteplength=1
> weightfactor=1
> for(i in 1:days)
> {
> #This is the stochastic component of the movement
> circularspace$stochasticxvel[day=i+1]=runif(1,min=-1,max=1)
> circularspace$stochasticyvel[day=i+1]=runif(1,min=-1,max=1)
> #This is calculating the movment speed
> circularspace$xvelocity[day=i+1]=weightfactor*(circularspace
> $xvelocity[day=i])+ (1-weightfactor)*circularspace
> $stochasticxvel[day=i+1]
See below. That looks all wrong. should just be [i] not
[day=i]>
> circularspace$yvelocity[day=i+1]=weightfactor*(circularspace
> $yvelocity[day=i])+ (1-weightfactor)*circularspace
> $stochasticyvel[day=i+1]
> #This is calculating the new position
> circularspace$xcoord[day=i+1]=circularspace$xcoord[day=i]
> +circularspace$xvelocity[day=i]/timesteplength
^^^^
here you need to learn to use <- for
assignment>
> circularspace$ycoord[day=i+1]=circularspace$ycoord[day=i]
> +circularspace$yvelocity[day=i]/timesteplength
> #Tis is checking the distance from the wall
> circularspace$xdistwall[day=i+1]=xmax-abs(circularspace$xcoord[day=i])
> circularspace$ydistwall[day=i+1]=ymax-abs(circularspace$ycoord[day=i])
I would have expected to see code that checked whether the radial
distance were greater than 10,
To do so would require calculating x^ + y^2. At the moment, you appear
to have a square bounding box.
And why is the distance on day 5 calculated in terms of day 4?
And you need to look at the definitions of logical operators. "=" is
not a logical operator. Above you were making unnecessary assignments,
now you are conflating "=" , one of the assignment operators, with
"==", the logical test for equality.
?"==" # or
?Comparison
>
> #This is updating the movement if distance to wall too small
> if (circularspace$xdistwall[day=i+1] <= mindist){
> circularspace$wallxvel[day=i+1]= -1*(circularspace$xcoord[day=i+1])
> circularspace$xvelocity[day=i+1]=weightfactor*circularspace
> $xvelocity[day=i]+ (1-weightfactor)*circularspace
> $stochasticxvel[day=i+1]+ circularspace$wallxvel[day=i+1]
> circularspace$xcoord[day=i+1]=circularspace$xcoord[day=i]
> +circularspace$xvelocity[day=i]/timesteplength
> }
> if (circularspace$ydistwall[day=i+1] <= mindist){
> circularspace$wallyvel[day=i+1]= -1*(circularspace$ycoord[day=i+1])
> circularspace$yvelocity[day=i+1]=weightfactor*circularspace
> $yvelocity[day=i]+ (1-weightfactor)*circularspace
> $stochasticyvel[day=i+1]+ circularspace$wallyvel[day=i+1]
> circularspace$ycoord[day=i+1]=circularspace$ycoord[day=i]
> +circularspace$yvelocity[day=i]/timesteplength
> }
> points(circularspace$xcoord[day=i+1],circularspace$ycoord[day=i+1])
> symbols(0,0,circles=radius,inches=F,add=T)
> symbols(0,0,circles=radius-2,inches=F,add=T)
> }
> circularspace
You might want to look at the worked example here:
http://bm2.genes.nig.ac.jp/RGM2/R_current/library/animation/man/brownian.motion.html
>
> Dr. Juliane Struve
> Environmental Scientist
> 10, Lynwood Crescent
> Sunningdale SL5 0BL
> 01344 620811
>
>
>
>
> ______________________________________________
> R-help at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.