I'm trying to get code to show the increase of density from the center outward.

I'm trying to get a code to show the increase of density from the center outward. The density increase 2kg/m^3 every meter outward from the center. The center is 3kg/M^3. I have working for a cone but don't know how to make it work for a sphere. I need it go from the center outward instead of the cone tip upward. Radius of the sphere is 5 meters. Code for the cone is below:

% Parameters for the cone
cone_height = 5; % in meters
cone_radius = 4; % in meters
% Density parameters
delta_0 = 3; % base density at the tip in kg/m^3
k = 2; % density increase rate in kg/m^3 per meter
% Grid of points
[z, theta] = meshgrid(linspace(0, cone_height, 50), linspace(0, 2*pi, 50));
% Radius at each height
r = (cone_radius / cone_height) * z;
% Polar coordinates to Cartesian coordinates
X = r .* cos(theta);
Y = r .* sin(theta);
% Density at each height
density = delta_0 + k * z;
% Create the figure
figure;
% Color gradient for density
surf(X, Y, z, density, 'FaceAlpha', 0.8);
colormap('jet'); 
colorbar; 
title('Shape 1 Visual');
xlabel('X (m)');
ylabel('Y (m)');
zlabel('Height (m)');
% Adjust view
view(-30, 30);
axis equal; 

[Image Analyst]

I'm not sure what it would look like. To have the density vary radially you'd need a volumetric (3-D) image. But then you'd only see the outer shell, which would all be at the same radius and thus the same color. Please mock up a picture of what you need to achieve, and give the "use case" (context).

[DGM]

While density may vary throughout a 3D solid, at least in the given example, the underlying relationship is a simple linear function of only one variable. It could be illustrated with a simple line plot, and it would be much easier to accurately read the plot.

Answer 1

[David Goodmanson]

Hi Matthew,

you could try a cutaway view (which doesn't represent the sphere as a solid)

n = 40;
r0 = 5;
[x y z] = sphere(n);
x1 = r0*x(n/2+1:end,:);
y1 = r0*y(n/2+1:end,:);
z1 = r0*z(n/2+1:end,:);
z2 = 0*z1;
r = sqrt(x1.^2+y1.^2+z1.^2);
rho = (3 + 2*r);
fig(1)
u = .4;
surf(x1,y1,z1,rho,'facealpha',1/2,'edgecolor',[u u u])
hold on
r = sqrt(x1.^2+y1.^2+z2.^2);
rho = (3 + 2*r);
surf(x1,y1,z2,rho,'facealpha',1/2,'edgecolor',[u u u])
hold off
axis equal
view([1 1 -1])
colorbar