we are on a good way, maybe...

This commit is contained in:
tTh
2022-12-05 13:10:40 +01:00
parent 8c9625b7df
commit 5b6df523fc
5 changed files with 104 additions and 43 deletions

View File

@@ -9,11 +9,12 @@ module realfield
implicit none
!-----------------------------------------------------------------------
!-
! definition of structures
!-
type massbody
real :: posx, posy
real :: heading = 0.21
real :: posx = 0, posy = 0
real :: heading = 0.29
real :: speed = 1.017
real :: mass = 1.0
integer :: serial = 666
@@ -22,63 +23,76 @@ end type
!-----------------------------------------------------------------------
contains
!-----------------------------------------------------------------------
subroutine barycentre_bodies(astres)
subroutine compute_barycentre_bodies(astres, bcx, bcy)
type(massbody), intent(in) :: astres(:)
real :: cx, cy
real, intent(out) :: bcx, bcy
integer :: foo
real :: cx, cy
!-
! May be we have to use DOUBLE RPECSION here ?
!-
cx = 0.0
cy = 0.0
do foo=1, ubound(astres, 1)
cx = cx + astres(foo)%posx
cy = cy + astres(foo)%posy
enddo
cx = cx / real(ubound(astres, 1))
cy = cy / real(ubound(astres, 1))
print *, 'barycentre:', cx, cy
bcx = cx / real(ubound(astres, 1))
bcy = cy / real(ubound(astres, 1))
end subroutine
!-----------------------------------------------------------------------
subroutine print_barycentre_bodies(astres)
type(massbody), intent(in) :: astres(:)
real :: cx, cy
call compute_barycentre_bodies(astres, cx, cy)
print *, "barycentre : ", cx, cy
end subroutine
!-----------------------------------------------------------------------
!-
! make a few solid body to play with...
! make a few solid body to play with...
!-
! planets : an array of type(massbody) to be filled
! coef : for setting the mass of the body
! sx, sy : borders of the universe
!-
subroutine create_some_planets(planets, coef, sx, sy)
type(massbody), intent(inout) :: planets(:)
real, intent(in) :: coef
integer, intent(in) :: sx, sy
integer :: foo
character(100) :: fmt
fmt = "(I4, ' | ', 2(F10.2, ' '), ' | ', 2F9.3, ' ', e12.3, I7)"
do foo=1, ubound(planets, 1)
!-
! the first planet is the home of Johnny Root
!-
if (foo .EQ. 1) then
!-
! the first planet is the home of Johnny Root
!-
planets(1)%posx = sx / 2
planets(1)%posy = sy / 2
planets(1)%mass = 37e8
planets(1)%mass = 29e8
planets(1)%serial = 1337
planets(1)%speed = 6.666
else
!-
! others are planets for peones
!-
planets(foo)%posx = rand() * real(sx-1)
planets(foo)%posy = rand() * real(sy-1)
planets(foo)%mass = 7e6 + coef*foo
planets(foo)%heading = 3.14159 * rand()
if (rand() .LT. 0.01) planets(foo)%speed = 2.718
if (rand() .LT. 0.08) planets(foo)%speed = 3.14159
planets(foo)%serial = foo*2 + 120
endif
write (*, fmt) foo, planets(foo)
enddo
end subroutine
!-----------------------------------------------------------------------
!-
! the basis of the kluge
!-
function compute_gravity(fx, fy, body)
real, intent(in) :: fx, fy
type(massbody), intent(in) :: body
@@ -88,7 +102,7 @@ function compute_gravity(fx, fy, body)
rx = fx - body%posx
ry = fy - body%posy
dist = sqrt( (rx*rx) + (ry*ry) )
if (dist .LT. 0.11) then
if (dist .LT. 0.08) then
! write (0, *) "dist too small ", dist
compute_gravity = 0e0
else