% Conductivity Excercise % Insert input parameteres (columns: Depth, porosity, b, conductivity; rows: % roks types (sand, siltstone, clay, marl,limestone) Tab=[1300,0.35,0.284,1.79;3500,0.30,0.379,1.58;4800,0.25,0.293,1.43;7200,0.15,0.298,1.68;10000,0.25,0.396,2.3]; % Define the water conductivity Kw=0.6; % Define the depth of the stratigraphy from 0 to 10 km z=[0:100:10000]; % Define the depth of the lithology change Depth=[0;Tab(:,1)]; % Initialize the matrix of the conductvity MagK1=zeros(100,1); MagK2=zeros(100,1); % Initialize the matrix of the thermal resistance mean MagThermResMean1=zeros(5,1); MagThermResMean2=zeros(5,1); % Initialize the matrix of the harmonic mean MagHarmM1=zeros(5,1); MagHarmM2=zeros(5,1); % Create a figure figure for i=1:5 NewDepth=(z(find(z==Depth(i)))+100:100:z(find(z==Depth(i+1))))./1000; % Estimate the porosity variation (function of depth) fi=Tab(i,2).*exp(-Tab(i,3).*NewDepth); % Estimate the thermal conductivity of each rock K1 (geometric mean) and K2 (square % root mean) K1=(Tab(i,4).^(1-fi)).*Kw.^fi; K2=((1-fi).*(Tab(i,4)).^0.5+fi.*Kw^0.5).^2; % plot the rocks' conductivities plot(K1,NewDepth.*-1,'Marker','diamond','LineWidth',1) hold on plot(K2,NewDepth.*-1,'Marker','v','LineWidth',1) hold on % Create xlabel xlabel('Thermal Conductivity (W/mK)','FontSize',20,'FontWeight','bold'); % Create ylabel ylabel('Depth (km)','FontSize',20,'FontWeight','bold'); % Estimate the aritmetic mean of the thermal conductivity for each % lithology K1Mean=mean(K1); K2Mean=mean(K2); % Estimate the fracion of each lithology and of the thermal resistance Fraction=(Depth(i+1)-Depth(i))./Depth(6); ThermResMean1=Fraction./K1Mean; ThermResMean2=Fraction./K2Mean; MagThermResMean1(i)=ThermResMean1; MagThermResMean2(i)=ThermResMean2; % Estimate the fracion of the harmonic mean HarmM1=((Depth(i+1)-Depth(i))/K1Mean); MagHarmM1(i)=HarmM1; HarmM2=((Depth(i+1)-Depth(i))/K2Mean); MagHarmM2(i)=HarmM2; % Write in two vectors the conductivities of the rocks Z=find(MagK1==0); S1=size(K1); MagK1(Z(1):Z(1)-1+S1(2))=K1; MagK2(Z(1):Z(1)-1+S1(2))=K2; end % Estimate the thermal resistance mean ThermResMean1=sum(MagThermResMean1); ThermResMean2=sum(MagThermResMean2); % Estimate the harmonic mean HarmonMean1=(Tab(5,1)./sum(MagHarmM1)); HarmonMean2=(Tab(5,1)./sum(MagHarmM2)); %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Ecercise 2 % Construct the geotherm with a variable gradient (from 30C/km to 21C/km at a depth of 10 km) % Initialize the matrix and associate the first two values T0=0; T1=30; MagT1=zeros(11,1); MagT1(1)=T0; MagT1(2)=T1; % Construct the geotherm every km for h=1:9 T1=T1+(30-h); MagT1(h+2)=T1; end % Interpolate the values of temperatures every 100 m (z) Depthkm=[0:1000:10000]; NewTemp=interp1(Depthkm,MagT1,z); % Estimate the thermal conductivity of the water for the two conditions if NewTemp<=137 Kw1=0.5648+1.878*10^-3.*NewTemp-7.231.*10^-6*(NewTemp).^2; else Kw1=0.6020+1.309*10^-3.*NewTemp-5.140.*10^-6*(NewTemp).^2; end % Estimate the conductivity of the matrix figure for j=1:5 % Estimate the thermal conductivity of the rocks'matrix as function of % temperature Km=1.8418+(Tab(j,4)-1.8418).*(1./(0.002732.*NewTemp+0.7463)-0.2485); % Define the depth of each rock (in m and km) with a step of 100 m DepthRock=[z(find(z==Depth(j)))+100:100:z(find(z==Depth(j+1)))]; NewDepthRock=DepthRock./1000; % Estimate the size of the vector of the depth of each rock SizeDR=size(DepthRock); SizeDR1=SizeDR(2); % Initialize the matrices of the thermal conductivity (with the size of the % vector of the depth of each rock) and associate the values of the thermal % conductivity of the rock matrix and water MagKm1=zeros(SizeDR1,1); MagKw1=zeros(SizeDR1,1); for l=1:101 PosDepth=find(z(l)==DepthRock); if PosDepth ~=0 Km1=Km(PosDepth); MagKm1(PosDepth)=Km1; Kw2=Kw1(PosDepth); MagKw1(PosDepth)=Kw2; else end % Estimate the porosity variation (function of depth) fi=Tab(j,2).*exp(-Tab(j,3).*NewDepthRock'); % Estimate the thermal conductivity of each rock K1 (geometric mean) and K2 (square % root mean) Kmw1=(MagKm1.^(1-fi)).*MagKw1.^fi; Kmw2=((1-fi).*(MagKm1).^0.5+fi.*MagKw1.^0.5).^2; end % plot the thermal conductivity (geometric mean and square root mean) plot(Kmw1,NewDepthRock.*-1,'Marker','<','LineWidth',1) hold on plot(Kmw2,NewDepthRock.*-1,'Marker','>','LineWidth',1) hold on % Create xlabel xlabel('Thermal Conductivity (W/mK)','FontSize',20,'FontWeight','bold'); % Create ylabel ylabel('Depth (km)','FontSize',20,'FontWeight','bold'); end