AfarTJ-22MyrTo10Myr.nb

Several points along the (mostly) non-extended edges of the Arabian, Nubian, and Somalian plates at the Afar triple junction have been selected and grouped with the rift that they are adjacent to:  Red Sea rift, East African rift, Aden rift.
These points are moved back to where they were 22 million years ago, or about the time the East Africa rift started spreading.
We assume that the current angular velocity vectors for the Arabian, Nubian, and Somalian plates have been constant for the past 22 Myr and will continue to be constant in the next few million years, during the time we are interested in modeling.

Revised April 23, 2019

Input

Select the future time for the model

modelTimeF=5;

Select the Euler pole data to be used for each of the three plates

poleLatSM=54.6131;
poleLongSM=-71.5876;
angSpeedSM=0.225471;

poleLatNB=50.3991;
poleLongNB=-36.5919;
angSpeedNB=0.189645;

poleLatAR=42.155;
poleLongAR=-0.270986;
angSpeedAR=0.456179;

Input the geographic data marking the edge of (mostly) undeformed continental crust for each of the three plates

The following matrices contain the geographic coordinates of several points along the (mostly) non-extended edges of the plate referenced in the matrix name.  Each record corresponds to a different boundary point.  The first column contains decimal latitudes (°) and the second column contains decimal longitudes (°).  The data were collected by visual interpretation of a Mercator map projection of a global DEM presented by GeoMapApp (www.geomapapp.org).

nubiaPtsToday=(({{12.662, 39.507}, {11.613, 39.590}, {10.367, 39.814}, {9.674, 39.744}, {8.69657, 39.4531}, {7.0467, 37.8161}}));

noRecsNB=6;

somaliaPtsToday=(({{7.2350, 38.8180}, {8.8330, 40.5044}, {9.440, 41.493}, {9.518, 42.047}, {9.748, 43.128}}));

noRecsSM=5;

arabiaPtsToday=(({{17.095, 43.124}, {15.441, 43.370}, {14.179, 43.576}, {13.693, 45.035}, {13.885, 45.756}, {14.218, 46.986}}));

noRecsAR=6;

User-Defined Functions

The function unitVect3D operates on a 3D vector (vect_) and outputs the Cartesian coordinates of the corresponding unit vector.

unitVect3D[vect_]:={(vect[[1]]/Norm[vect]),(vect[[2]]/Norm[vect]),(vect[[3]]/Norm[vect])};

zRotation[w_,dT_]:={{Cos[(w dT)Degree],-Sin[(w dT)Degree],0},{Sin[(w dT)Degree],Cos[(w dT)Degree],0},{0,0,1}};

The function geog2Cart operates on two geographic coordinates (lat_ and long_) and outputs the Cartesian coordinates of the unit location vector to that point.

geog2Cart[lat_, long_] := {Cos[lat Degree] Cos[long Degree], Cos[lat Degree] Sin[long Degree], Sin[lat Degree]};

The function cart2Geog operates on a 3D unit location vector (vect_) and outputs the geographic coordinates of that point in latitude and longitude.

cart2Geog[vect_]:=Module[{lat,long,a,b,c,d,e,f},a=ArcSin[vect[[3]]];b={vect[[1]],vect[[2]],0};c=If[((Abs[vect[[1]]]<(1×(10^(-14))))&&(Abs[vect[[2]]]<(1×(10^(-14))))),{1,1,0},{vect[[1]]/Norm[b],vect[[2]]/Norm[b],0}];d={1,0,0};e=VectorAngle[c,d];f=If[(vect[[2]]<0),(-e),(e)];lat=a(180/π);long=If[((Abs[vect[[1]]]<(1×(10^(-14))))&&(Abs[vect[[2]]]<(1×(10^(-14))))),0,(f(180/π))];{lat,long}];

The function distAzFromVectors operates on a 3D unit location vector to a reference point (refPtVect_) and a 3D unit location vector to another point (otherPtVect_), and outputs the distance in km and azimuth in degrees from the reference point to the other point.

distAzFromVectors[refPtVect_,otherPtVect_]:=Module[{radEarth,circEarth,northPole,θ1,θ2,θ3,θ4,normVect1,normVect2,normVect3,distKm,bearing,answer},radEarth=6371.02;circEarth=2*π*radEarth;northPole={0,0,1};normVect1 = Cross[northPole,refPtVect];
normVect2 = Cross[otherPtVect,refPtVect];θ1=VectorAngle[refPtVect,otherPtVect];
distKm=(θ1/(2*π))*circEarth;θ2=(180/π)*(VectorAngle[normVect1,normVect2]);θ3=(180/π)*(VectorAngle[normVect1,otherPtVect]);normVect3=Cross[refPtVect,normVect1];θ4=(180/π)*(VectorAngle[normVect3,otherPtVect]);bearing=If[((θ3==90)∨(θ3==180)),If[(θ4<=90),0,180],If[(θ3>90),(360-θ2),θ2]];answer={distKm,bearing};answer];

The function circMotion operates on a 3D unit location vector (x_), a coordinate transformation matrix (m1_), an angular speed (w_) and a time interval (dT_) in Myr, and outputs the unit location vector to that point after it rotates around its axis.

circMotion[x_,m1_,w_,dT_]:=Module[{m2,m3,answer},m2=zRotation[w,dT];m3=Inverse[m1];answer=m3.m2.m1.x;answer];

The function spherePlot1 operates on a set of 3D unit location vectors of computed data (inData_) and two sets of 3D unit location vectors of colored markers (markers1_ and markers2_) to help establish the coordinate system(s), and outputs a graphic of a unit sphere on which the input points are plotted.

spherePlot1[inData_,markers1_,markers2_]:=Module[{g1,g2,g3,g4,answer},g1=Graphics3D[Sphere[{0,0,0},1],AspectRatio->1,BoxRatios->{1,1,1},PlotRange->All,PlotRangePadding->0.1,ColorOutput->GrayLevel,Lighting->"Neutral"];g2=ListPointPlot3D[inData,AspectRatio->1,BoxRatios->{1,1,1},PlotStyle->Blue,PlotRange->All,PlotRangePadding->0.1];g3=ListPointPlot3D[markers1,AspectRatio->1,BoxRatios->{1,1,1},PlotStyle->Red,PlotRange->All,PlotRangePadding->0.1];g4=ListPointPlot3D[markers2,AspectRatio->1,BoxRatios->{1,1,1},PlotStyle->Green,PlotRange->All,PlotRangePadding->0.1];answer=Show[g1,g2,g3,g4];answer];

cycloid[obsΩ_,movΩ_,refPt_,t_]:=Module[{a,b,c,d,e,f,g,h,i,j,k,l,m1, m2, m3, m4, m5,n,p,answer},a={1,0,0};b={0,1,0};c={0,0,1};d=unitVect3D[obsΩ];e=Norm[obsΩ];f=unitVect3D[movΩ];g=Norm[movΩ];h=unitVect3D[d\[Cross]f];i=h;j=unitVect3D[h\[Cross]d];k=unitVect3D[i\[Cross]f];m1={{k[[1]],k[[2]],k[[3]]},{i[[1]],i[[2]],i[[3]]}, {f[[1]],f[[2]],f[[3]]}};m2={{Cos[(-g Degree) t],Sin[(-g Degree) t],0},{-Sin[(-g Degree) t],Cos[(-g Degree) t],0},{0,0,1}};l=VectorAngle[d,k];If[l<(π/2),(n=1),(n=(-1))];p = (n) (VectorAngle[f,d]);m3={{Cos[p],0,-Sin[p]},{0,1,0},{Sin[p],0,Cos[p]}};m4={{Cos[(e Degree) t],Sin[(e Degree) t],0},{-Sin[(e Degree) t],Cos[(e Degree) t],0},{0,0,1}};m5={{j[[1]],h[[1]],d[[1]]},{j[[2]],h[[2]],d[[2]]},{j[[3]],h[[3]],d[[3]]}};answer=m5.m4.m3.m2.m1.refPt;answer];

The user-defined function northArrow takes the average of all the unit location vectors of points along the plate boundary (viewPoint) and a user-defined scale (nArrowScale), generally ranging from about 0.05 to 0.1, and finds a suitable endpoint for a vector pointing toward north in the output graphics.

northArrow[viewPoint_,nArrowScale_]:=Module[{north,center,normal,northward,answer},north={0,0,1};center={0,0,0};normal=Cross[north,viewPoint];northward=Cross[viewPoint,normal]*nArrowScale;answer={center,northward};answer];

midpoint[cartVect1_,cartVect2_]:=Module[{a,answer},a=(cartVect1+cartVect2)/2;answer={(a[[1]]/Norm[a]),(a[[2]]/Norm[a]),(a[[3]]/Norm[a])};answer];

intersection[vect1_,vect2_,vect3_,vect4_]:=Module[{a,b,c,d,answer},a=Cross[vect1,vect2];b=Cross[vect3,vect4];c=Cross[b,a];
d={(c[[1]]/Norm[c]),(c[[2]]/Norm[c]),(c[[3]]/Norm[c])};answer=If[(VectorAngle[d,vect2]>(π/2)),-d,d];answer];

average3Vectors[vect1_,vect2_,vect3_]:=Module[{a,b,answer},a=(vect1+vect2+vect3)/3;answer={(a[[1]]/Norm[a]),(a[[2]]/Norm[a]),(a[[3]]/Norm[a])};answer];

addExplanation[plateName_,geogCoords_,modTime_]:=Module[{a,answer},a=Table[Prepend[geogCoords[[i]],plateName],{i,1,Length[geogCoords]}];answer=Table[Append[a[[i]],modTime],{i,1,Length[a]}];answer];

Computation

The model time is 22 Myr ago (-22).

modelTime=-22;

Move each set of points back to where they would have been at the selected time with a circular motion around the respective plate pole.

Define the axes of the geographic coordinate system

xGeog={1,0,0};

yGeog={0,1,0};

zGeog={0,0,1};

Nubian Plate

Convert the input pole location data for the pole to the Nubian plate to Cartesian coordinates

poleEulerNB=geog2Cart[poleLatNB,poleLongNB];

Define the angular speed “ω” for the Nubian plate

ωNB=angSpeedNB;

Define the axes of the Euler pole coordinate system

zPoleNB=poleEulerNB;

xPoleNB=Cross[zPoleNB,zGeog];

yPoleNB=Cross[zPoleNB,xPoleNB];

Find the Cartesian coordinates of the reference points on the Somalian plate

cartNBtoday=Table[geog2Cart[nubiaPtsToday[[i,1]],nubiaPtsToday[[i,2]]],{i,1,noRecsNB}];

Define the transformation matrix between the geographic coordinate system and the Euler pole coordinate system

jNB=(({{Dot[xGeog,xPoleNB], Dot[yGeog,xPoleNB], Dot[zGeog,xPoleNB]}, {Dot[xGeog,yPoleNB], Dot[yGeog,yPoleNB], Dot[zGeog,yPoleNB]}, {Dot[xGeog,zPoleNB], Dot[yGeog,zPoleNB], Dot[zGeog,zPoleNB]}}));

Rotate the reference point around the Euler pole by an angle equal to the angular speed multiplied by the time interval.

rotPtsNB22mybp=Table[circMotion[cartNBtoday[[i]],jNB,ωNB,modelTime],{i,1,noRecsNB}];

rotGeogNB22mybp=Table[cart2Geog[rotPtsNB22mybp[[i]]],{i,1,Length[rotPtsNB22mybp]}];

MatrixForm[%]

Somalian Plate

Convert the input pole location data for the pole to the Somalian plate to Cartesian coordinates

poleEulerSM=geog2Cart[poleLatSM,poleLongSM];

Define the angular speed “ω” for the Somalian plate

ωSM=angSpeedSM;

Define the axes of the Euler pole coordinate system

zPoleSM=poleEulerSM;

xPoleSM=Cross[zPoleSM,zGeog];

yPoleSM=Cross[zPoleSM,xPoleSM];

Find the Cartesian coordinates of the reference points on the Somalian plate

cartSMtoday=Table[geog2Cart[somaliaPtsToday[[i,1]],somaliaPtsToday[[i,2]]],{i,1,noRecsSM}];

Define the transformation matrix between the geographic coordinate system and the Euler pole coordinate system

jSM=(({{Dot[xGeog,xPoleSM], Dot[yGeog,xPoleSM], Dot[zGeog,xPoleSM]}, {Dot[xGeog,yPoleSM], Dot[yGeog,yPoleSM], Dot[zGeog,yPoleSM]}, {Dot[xGeog,zPoleSM], Dot[yGeog,zPoleSM], Dot[zGeog,zPoleSM]}}));

Rotate the reference point around the Euler pole by an angle equal to the angular speed multiplied by the time interval.

rotPtsSM22mybp=Table[circMotion[cartSMtoday[[i]],jSM,ωSM,modelTime],{i,1,noRecsSM}];

rotGeogSM22mybp=Table[cart2Geog[rotPtsSM22mybp[[i]]],{i,1,Length[rotPtsSM22mybp]}];

MatrixForm[%]

Arabian Plate

Convert the input pole location data for the pole to the Arabian plate to Cartesian coordinates

poleEulerAR=geog2Cart[poleLatAR,poleLongAR];

Define the angular speed “ω” for the Arabian plate

ωAR=angSpeedAR;

Define the axes of the Euler pole coordinate system

zPoleAR=poleEulerAR;

xPoleAR=Cross[zPoleAR,zGeog];

yPoleAR=Cross[zPoleAR,xPoleAR];

Find the Cartesian coordinates of the reference points on the Arabian plate

cartARtoday=Table[geog2Cart[arabiaPtsToday[[i,1]],arabiaPtsToday[[i,2]]],{i,1,noRecsAR}];

Define the transformation matrix between the geographic coordinate system and the Euler pole coordinate system

jAR=(({{Dot[xGeog,xPoleAR], Dot[yGeog,xPoleAR], Dot[zGeog,xPoleAR]}, {Dot[xGeog,yPoleAR], Dot[yGeog,yPoleAR], Dot[zGeog,yPoleAR]}, {Dot[xGeog,zPoleAR], Dot[yGeog,zPoleAR], Dot[zGeog,zPoleAR]}}));

Rotate the reference point around the Euler pole by an angle equal to the angular speed multiplied by the time interval.

rotPtsAR22mybp=Table[circMotion[cartARtoday[[i]],jAR,ωAR,modelTime],{i,1,noRecsAR}];

rotGeogAR22mybp=Table[cart2Geog[rotPtsAR22mybp[[i]]],{i,1,Length[rotPtsAR22mybp]}];

MatrixForm[%]

Construct tables of the current and back-rotated points along the continental margins of the Nubian, Somalian, and Arabian Plates

modTime=0;

plateName="NubiaContEdge";

nubiaContEdgeToday=addExplanation[plateName,nubiaPtsToday,modTime];

plateName="SomaliaContEdge";

somaliaContEdgeToday=addExplanation[plateName,somaliaPtsToday,modTime];

plateName="ArabiaContEdge";

arabiaContEdgeToday=addExplanation[plateName,arabiaPtsToday,modTime];

Clear[plateName];

resultsContEdgeToday=Join[nubiaContEdgeToday,somaliaContEdgeToday,arabiaContEdgeToday];

modTime=-22;

plateName="NubiaContEdge";

nubiaContEdge22mybp=addExplanation[plateName,rotGeogNB22mybp,modTime];

plateName="SomaliaContEdge";

somaliaContEdge22mybp=addExplanation[plateName,rotGeogSM22mybp,modTime];

plateName="ArabiaContEdge";

arabiaContEdge22mybp=addExplanation[plateName,rotGeogAR22mybp,modTime];

Clear[plateName];

resultsContEdge22mybp=Join[nubiaContEdge22mybp,somaliaContEdge22mybp,arabiaContEdge22mybp];

Rotate the current continental margins forward

rotPtsNB5Myr=Table[circMotion[cartNBtoday[[i]],jNB,ωNB,modelTimeF],{i,1,noRecsNB}];

rotGeogNB5Myr=Table[cart2Geog[rotPtsNB5Myr[[i]]],{i,1,Length[rotPtsNB5Myr]}];

rotPtsSM5Myr=Table[circMotion[cartSMtoday[[i]],jSM,ωSM,modelTimeF],{i,1,noRecsSM}];

rotGeogSM5Myr=Table[cart2Geog[rotPtsSM5Myr[[i]]],{i,1,Length[rotPtsSM5Myr]}];

rotPtsAR5Myr=Table[circMotion[cartARtoday[[i]],jAR,ωAR,modelTimeF],{i,1,noRecsAR}];

rotGeogAR5Myr=Table[cart2Geog[rotPtsAR5Myr[[i]]],{i,1,Length[rotPtsAR5Myr]}];

Construct table of the forward-rotated points along the continental margins of the Nubian, Somalian, and Arabian Plates

plateName="NubiaContEdge";

nubiaContEdge5Myr=addExplanation[plateName,rotGeogNB5Myr,modelTimeF];

plateName="SomaliaContEdge";

somaliaContEdge5Myr=addExplanation[plateName,rotGeogSM5Myr,modelTimeF];

plateName="ArabiaContEdge";

arabiaContEdge5Myr=addExplanation[plateName,rotGeogAR5Myr,modelTimeF];

Clear[plateName];

resultsContEdge5Myr=Join[nubiaContEdge5Myr,somaliaContEdge5Myr,arabiaContEdge5Myr];

Finding the Initial Triple Junction

Red Sea Rift

midpt1at22mybp=midpoint[rotPtsNB22mybp[[1]],rotPtsAR22mybp[[1]]];

midpt2at22mybp=midpoint[rotPtsNB22mybp[[3]],rotPtsAR22mybp[[3]]];

East Africa Rift

midpt3at22mybp=midpoint[rotPtsNB22mybp[[5]],rotPtsSM22mybp[[2]]];

midpt4at22mybp=midpoint[rotPtsNB22mybp[[6]],rotPtsSM22mybp[[1]]];

Aden Rift

midpt5at22mybp=midpoint[rotPtsSM22mybp[[5]],rotPtsAR22mybp[[6]]];

midpt6at22mybp=midpoint[rotPtsSM22mybp[[3]],rotPtsAR22mybp[[4]]];

Intersections and mean triple junction

The method of finding the triple junction use here is to find the mean of the three riftline intersections.  Consequently, no microplate is modeled to grow in this triple junction over time.

tj1=intersection[midpt4at22mybp,midpt3at22mybp,midpt2at22mybp,midpt1at22mybp];

tj2=intersection[midpt5at22mybp,midpt6at22mybp,midpt3at22mybp,midpt4at22mybp];

tj3=intersection[midpt1at22mybp,midpt2at22mybp,midpt6at22mybp,midpt5at22mybp];

meanTJcart22mybp=average3Vectors[tj1,tj2,tj3];

riftTJpoints22mybp={midpt1at22mybp,midpt2at22mybp,meanTJcart22mybp,midpt3at22mybp,midpt4at22mybp,midpt3at22mybp,meanTJcart22mybp,midpt6at22mybp,midpt5at22mybp};

nubiaPtsCart22mybp={midpt1at22mybp,midpt2at22mybp,meanTJcart22mybp,midpt3at22mybp,midpt4at22mybp};

nubiaPtsGeog22mybp=Table[cart2Geog[nubiaPtsCart22mybp[[i]]],{i,1,Length[nubiaPtsCart22mybp]}];

somaliaPtsCart22mybp={midpt4at22mybp,midpt3at22mybp,meanTJcart22mybp,midpt6at22mybp,midpt5at22mybp};

somaliaPtsGeog22mybp=Table[cart2Geog[somaliaPtsCart22mybp[[i]]],{i,1,Length[somaliaPtsCart22mybp]}];

arabiaPtsCart22mybp={midpt1at22mybp,midpt2at22mybp,meanTJcart22mybp,midpt6at22mybp,midpt5at22mybp};

arabiaPtsGeog22mybp=Table[cart2Geog[arabiaPtsCart22mybp[[i]]],{i,1,Length[arabiaPtsCart22mybp]}];

Construct table of model points at 22 Myr BP

modTime=-22;

plateName="nubia";

nubiaGeog22mybp=addExplanation[plateName,nubiaPtsGeog22mybp,modTime];

plateName="somalia";

somaliaGeog22mybp=addExplanation[plateName,somaliaPtsGeog22mybp,modTime];

plateName="arabia";

arabiaGeog22mybp=addExplanation[plateName,arabiaPtsGeog22mybp,modTime];

Clear[plateName];

resultsGeog22mybp=Join[nubiaGeog22mybp,somaliaGeog22mybp,arabiaGeog22mybp];

Prepare the graphic output

Parameters that affect the graphic output

Set the size of the north arrow on the graphic.

northArrowScale=0.1;

Choose the distance from Earth' s center (expressed in Earth radii) from  which the boundary will be viewed in the graphic.  The value assigned to this variable must be greater than 1.

viewDistance=1.28;

centerOfMap=geog2Cart[11.5,42.0];

Variable averageVV is the point from which the boundary is viewed looking toward the center of Earth, as plotted using the built-in Mathematica function Graphics3D

averageVV=viewDistance*centerOfMap;

Code that creates the initial graphic output

Note :  Some useful (easily visible) built - in colors for use in the graphics include the following :  Red, Green, Blue, Black, Gray, Cyan, Magenta, Brown, Orange, Pink, and Purple.

boundaryOut22mybp = Graphics3D[{Arrow[northArrow[averageVV,northArrowScale]],Black,Line[Table[rotPtsAR22mybp[[i]],{i,1,Length[rotPtsAR22mybp]}]],Line[Table[rotPtsSM22mybp[[i]],{i,1,Length[rotPtsSM22mybp]}]],Line[Table[rotPtsNB22mybp[[i]],{i,1,Length[rotPtsNB22mybp]}]],Red,Line[Table[riftTJpoints22mybp[[i]],{i,1,Length[riftTJpoints22mybp]}]],Blue, Dashed,Line[Table[cartARtoday[[i]],{i,1,Length[cartARtoday]}]], Line[Table[cartSMtoday[[i]],{i,1,Length[cartSMtoday]}]], Line[Table[cartNBtoday[[i]],{i,1,Length[cartNBtoday]}]]}, AspectRatio->Automatic, Boxed->False,ViewVector->{averageVV,{0,0,0}}];

boundaryOut22mybp

The graphic shown above shows the Afar triple junction (red) at -22 MYBP. 
• The original triple junction at 22 MYBP is shown by the solid red lines in the lower left, and the solid black lines on either side are the edges of the current non-extended continental crust as rotated back over 22 Myr.
• The blue dashed lines are the edges of the current non-extended continental crust today.

ClearAll[tj1,tj2,tj3];

Growing the triple junction from 22 MYBP to today

modelTime0=22;

Nubian Plate

The following table contains the latitude (first column) and longitude (second column) for the 5 model points on the Nubian plate 22 million years ago.  The third record is the modeled location of the triple junction at that time.

MatrixForm[nubiaPtsGeog22mybp]

Rotate the reference point around the Euler pole by an angle equal to the angular speed multiplied by the time interval.

rotPtsNB0=Table[circMotion[nubiaPtsCart22mybp[[i]],jNB,ωNB,modelTime0],{i,1,Length[nubiaPtsCart22mybp]}];

Somalian Plate

The following table contains the latitude (first column) and longitude (second column) for the 5 model points on the Somalian plate 22 million years ago.  The third record is the modeled location of the triple junction at that time.

MatrixForm[somaliaPtsGeog22mybp]

Rotate the reference point around the Euler pole by an angle equal to the angular speed multiplied by the time interval.

rotPtsSM0=Table[circMotion[somaliaPtsCart22mybp[[i]],jSM,ωSM,modelTime0],{i,1,Length[somaliaPtsCart22mybp]}];

Arabian Plate

The following table contains the latitude (first column) and longitude (second column) for the 5 model points on the Arabian plate 22 million years ago.  The third record is the modeled location of the triple junction at that time.

MatrixForm[arabiaPtsGeog22mybp]

Rotate the reference point around the Euler pole by an angle equal to the angular speed multiplied by the time interval.

rotPtsAR0=Table[circMotion[arabiaPtsCart22mybp[[i]],jAR,ωAR,modelTime0],{i,1,Length[arabiaPtsCart22mybp]}];

Red Sea Rift

midpt1Today=midpoint[rotPtsNB0[[1]],rotPtsAR0[[1]]];

midpt2Today=midpoint[rotPtsNB0[[2]],rotPtsAR0[[2]]];

East Africa Rift

midpt3Today=midpoint[rotPtsNB0[[4]],rotPtsSM0[[2]]];

midpt4Today=midpoint[rotPtsNB0[[5]],rotPtsSM0[[1]]];

Aden Rift

midpt5Today=midpoint[rotPtsSM0[[5]],rotPtsAR0[[5]]];

midpt6Today=midpoint[rotPtsSM0[[4]],rotPtsAR0[[4]]];

Intersections and mean triple junction

The method of finding the triple junction use here is to find the mean of the three riftline intersections.  Consequently, no microplate is modeled to grow in this triple junction over time.

tj1=intersection[midpt4Today,midpt3Today,midpt2Today,midpt1Today];

tj2=intersection[midpt5Today,midpt6Today,midpt3Today,midpt4Today];

tj3=intersection[midpt1Today,midpt2Today,midpt6Today,midpt5Today];

meanTJcartToday=average3Vectors[tj1,tj2,tj3];

riftTJpoints0={midpt1Today,midpt2Today,meanTJcartToday,midpt3Today,midpt4Today,midpt3Today,meanTJcartToday,midpt6Today,midpt5Today};

nubiaPtsCartToday={midpt1Today,midpt2Today,meanTJcartToday,midpt3Today,midpt4Today};

nubiaPtsGeogToday=Table[cart2Geog[nubiaPtsCartToday[[i]]],{i,1,Length[nubiaPtsCartToday]}];

The following table contains the latitude (first column) and longitude (second column) for the 5 model points on the Nubian plate today.  The third record is the modeled location of the triple junction today.

MatrixForm[nubiaPtsGeogToday]

somaliaPtsCartToday={midpt4Today,midpt3Today,meanTJcartToday,midpt6Today,midpt5Today};

somaliaPtsGeogToday=Table[cart2Geog[somaliaPtsCartToday[[i]]],{i,1,Length[somaliaPtsCartToday]}];

The following table contains the latitude (first column) and longitude (second column) for the 5 model points on the Somalian plate today.  The third record is the modeled location of the triple junction today.

MatrixForm[somaliaPtsGeogToday]

arabiaPtsCartToday={midpt1Today,midpt2Today,meanTJcartToday,midpt6Today,midpt5Today};

arabiaPtsGeogToday=Table[cart2Geog[arabiaPtsCartToday[[i]]],{i,1,Length[arabiaPtsCartToday]}];

The following table contains the latitude (first column) and longitude (second column) for the 5 model points on the Arabian plate today.  The third record is the modeled location of the triple junction today.

MatrixForm[arabiaPtsGeogToday]

Construct table of model points at their current locations

modTime=0;

plateName="nubia";

nubiaGeogToday=addExplanation[plateName,nubiaPtsGeogToday,modTime];

plateName="somalia";

somaliaGeogToday=addExplanation[plateName,somaliaPtsGeogToday,modTime];

plateName="arabia";

arabiaGeogToday=addExplanation[plateName,arabiaPtsGeogToday,modTime];

Clear[plateName];

resultsGeogToday=Join[nubiaGeogToday,somaliaGeogToday,arabiaGeogToday];

Prepare the graphic output

Code that creates the graphic output

Note :  Some useful (easily visible) built - in colors for use in the graphics include the following :  Red, Green, Blue, Black, Gray, Cyan, Magenta, Brown, Orange, Pink, and Purple.

boundaryOut22mybp = Graphics3D[{Arrow[northArrow[averageVV,northArrowScale]],Black,Line[Table[rotPtsAR0[[i]],{i,1,Length[rotPtsAR0]}]],Line[Table[rotPtsSM0[[i]],{i,1,Length[rotPtsSM0]}]],Line[Table[rotPtsNB0[[i]],{i,1,Length[rotPtsNB0]}]],Red,Line[Table[riftTJpoints0[[i]],{i,1,Length[riftTJpoints0]}]],Blue, Dashed,Line[Table[cartARtoday[[i]],{i,1,Length[cartARtoday]}]], Line[Table[cartSMtoday[[i]],{i,1,Length[cartSMtoday]}]], Line[Table[cartNBtoday[[i]],{i,1,Length[cartNBtoday]}]], Black,Dashed,Line[Table[rotPtsAR22mybp[[i]],{i,1,Length[rotPtsAR22mybp]}]],Line[Table[rotPtsSM22mybp[[i]],{i,1,Length[rotPtsSM22mybp]}]],Line[Table[rotPtsNB22mybp[[i]],{i,1,Length[rotPtsNB22mybp]}]],Red,Dashed,Line[Table[riftTJpoints22mybp[[i]],{i,1,Length[riftTJpoints22mybp]}]]}, AspectRatio->Automatic, Boxed->False,ViewVector->{averageVV,{0,0,0}}];

boundaryOut22mybp

The graphic shown above shows the Afar triple junction (red) at -22 MYBP and today. 
• The original triple junction at 22 MYBP is the dashed red line in the lower left, and the dashed black lines on either side are the edges of the current non-extended continental crust as rotated back over 22 Myr.
• The solid red lines in the middle of the diagram, with the black north arrow sticking out of them, is the current (model) triple junction.They are flanked by solid black lines, which are the displaced plate boundary from 22 MYBP (the red dashed lines).  The blue dashed lines are the edges of the current non - extended continental crust today.  The area between the solid black lines here is the new oceanic crust + extended continental crust of the last 22 Myr.

ClearAll[tj1,tj2,tj3];

Growing the triple junction from Today to a Future Time

modelTimeF

Nubian Plate

Rotate the reference point around the Euler pole by an angle equal to the angular speed multiplied by the time interval.

rotPtsNBF=Table[circMotion[nubiaPtsCartToday[[i]],jNB,ωNB,modelTimeF],{i,1,Length[nubiaPtsCartToday]}];

Somalian Plate

Rotate the reference point around the Euler pole by an angle equal to the angular speed multiplied by the time interval.

rotPtsSMF=Table[circMotion[somaliaPtsCartToday[[i]],jSM,ωSM,modelTimeF],{i,1,Length[somaliaPtsCartToday]}];

Arabian Plate

Rotate the reference point around the Euler pole by an angle equal to the angular speed multiplied by the time interval.

rotPtsARF=Table[circMotion[arabiaPtsCartToday[[i]],jAR,ωAR,modelTimeF],{i,1,Length[arabiaPtsCartToday]}];

Red Sea Rift

midpt1F=midpoint[rotPtsNBF[[1]],rotPtsARF[[1]]];

midpt2F=midpoint[rotPtsNBF[[2]],rotPtsARF[[2]]];

East Africa Rift

midpt3F=midpoint[rotPtsNBF[[4]],rotPtsSMF[[2]]];

midpt4F=midpoint[rotPtsNBF[[5]],rotPtsSMF[[1]]];

Aden Rift

midpt5F=midpoint[rotPtsSMF[[5]],rotPtsARF[[5]]];

midpt6F=midpoint[rotPtsSMF[[4]],rotPtsARF[[4]]];

Intersections and mean triple junction

The method of finding the triple junction use here is to find the mean of the three riftline intersections.  Consequently, no microplate is modeled to grow in this triple junction over time.

tj1=intersection[midpt4F,midpt3F,midpt2F,midpt1F];

tj2=intersection[midpt5F,midpt6F,midpt3F,midpt4F];

tj3=intersection[midpt1F,midpt2F,midpt6F,midpt5F];

meanTJcartF=average3Vectors[tj1,tj2,tj3];

riftTJpointsF={midpt1F,midpt2F,meanTJcartF,midpt3F,midpt4F,midpt3F,meanTJcartF,midpt6F,midpt5F};

nubiaPtsCartF={midpt1F,midpt2F,meanTJcartF,midpt3F,midpt4F};

nubiaPtsGeogF=Table[cart2Geog[nubiaPtsCartF[[i]]],{i,1,Length[nubiaPtsCartF]}];

The following table contains the latitude (first column) and longitude (second column) for the 5 model points on the Nubian plate after "modelTimeF" million years in the future.  The third record is the modeled location of the triple junction at that time.

MatrixForm[%]

somaliaPtsCartF={midpt4F,midpt3F,meanTJcartF,midpt6F,midpt5F};

somaliaPtsGeogF=Table[cart2Geog[somaliaPtsCartF[[i]]],{i,1,Length[somaliaPtsCartF]}];

The following table contains the latitude (first column) and longitude (second column) for the 5 model points on the Somalian plate after "modelTimeF" million years in the future.  The third record is the modeled location of the triple junction at that time.

MatrixForm[%]

arabiaPtsCartF={midpt1F,midpt2F,meanTJcartF,midpt6F,midpt5F};

arabiaPtsGeogF=Table[cart2Geog[arabiaPtsCartF[[i]]],{i,1,Length[arabiaPtsCartF]}];

The following table contains the latitude (first column) and longitude (second column) for the 5 model points on the Arabian plate after "modelTimeF" million years in the future.  The third record is the modeled location of the triple junction at that time.

MatrixForm[%]

Construct table of model points at a specified future time

modTime=5;

plateName="nubia";

nubiaGeogF=addExplanation[plateName,nubiaPtsGeogF,modTime];

plateName="somalia";

somaliaGeogF=addExplanation[plateName,somaliaPtsGeogF,modTime];

plateName="arabia";

arabiaGeogF=addExplanation[plateName,arabiaPtsGeogF,modTime];

Clear[plateName];

resultsGeogF=Join[nubiaGeogF,somaliaGeogF,arabiaGeogF];

Prepare the graphic output

Code that creates the graphic output

Note :  Some useful (easily visible) built - in colors for use in the graphics include the following :  Red, Green, Blue, Black, Gray, Cyan, Magenta, Brown, Orange, Pink, and Purple.

boundaryOutF = Graphics3D[{Arrow[northArrow[averageVV,northArrowScale]],Black,Line[Table[rotPtsARF[[i]],{i,1,Length[rotPtsARF]}]],Line[Table[rotPtsSMF[[i]],{i,1,Length[rotPtsSMF]}]],Line[Table[rotPtsNBF[[i]],{i,1,Length[rotPtsNBF]}]],Red,Line[Table[riftTJpointsF[[i]],{i,1,Length[riftTJpointsF]}]],Blue, Dashed,Line[Table[cartARtoday[[i]],{i,1,Length[cartARtoday]}]], Line[Table[cartSMtoday[[i]],{i,1,Length[cartSMtoday]}]], Line[Table[cartNBtoday[[i]],{i,1,Length[cartNBtoday]}]], Black,Dashed,Line[Table[rotPtsAR0[[i]],{i,1,Length[rotPtsAR0]}]],Line[Table[rotPtsSM0[[i]],{i,1,Length[rotPtsSM0]}]],Line[Table[rotPtsNB0[[i]],{i,1,Length[rotPtsNB0]}]],Line[Table[rotPtsAR22mybp[[i]],{i,1,Length[rotPtsAR22mybp]}]],Line[Table[rotPtsSM22mybp[[i]],{i,1,Length[rotPtsSM22mybp]}]],Line[Table[rotPtsNB22mybp[[i]],{i,1,Length[rotPtsNB22mybp]}]],Red,Dashed,Line[Table[riftTJpoints0[[i]],{i,1,Length[riftTJpoints0]}]],Line[Table[riftTJpoints22mybp[[i]],{i,1,Length[riftTJpoints22mybp]}]]}, AspectRatio->Automatic, Boxed->False,ViewVector->{averageVV,{0,0,0}}];

boundaryOutF

The graphic shown above shows the Afar triple junction (red) at -22 MYBP, today, and 10 Myr in the future.  
• The original triple junction at 22 MYBP is the dashed red line in the lower left, and the dashed black lines on either side are the edges of the current non-extended continental crust as rotated back over 22 Myr.
• The red dashed lines in the middle of the diagram, with the black north arrow sticking out of them, is the current (model) triple junction.They are flanked by black dashed lines, which are the displaced plate boundary from 22 MYBP.  The blue dashed lines are the edges of the current non - extended continental crust today.  The area between the dashed black lines here is the new oceanic crust + extended continental crust of the last 22 Myr.
• The solid red lines to the upper right is the ridge - triple - junction system 10 Myr in the future.The solid black lines are the displaced ridge - triple - junction system of today as displaced on each of the three plates over 10 Myr of finite displacement.

ClearAll[tj1,tj2,tj3];

Construct the consolidated table of model points

outFile1=Join[resultsGeog22mybp,resultsGeogToday,resultsGeogF,resultsContEdge22mybp,resultsContEdgeToday,resultsContEdge5Myr];

MatrixForm[%]

Export["/Users/vince/Desktop/Ryley/PlotData4GeoMapApp.xlsx",outFile1]