Code to check whether an input is a repfigit (Keith number).

f[list_]:=Append[Rest[list],Apply[Plus,list]]
knumsums[list_]:=With[{val=FromDigits[list]},Take[NestWhile[f,list,Last[#]<val&],-2]]
knum[list_]:=Last[knumsums[list]]===FromDigits[list]

Code to generate the family of repfigit equations for a given number of digits. Solutions to any of these equations give Keith numbers of that size. A given equation might have zero, one, or more than one solution over the integers.

keithEquations[len_Integer/;len>0]:=Module[{matrow,n,list,res,vecs},res=list[];
Do[matrow[j]=Table[KroneckerDelta[k,j+1],{k,len}],{j,len-1}];
matrow[len]=Table[1,{len}];
n=len;
While[9*Apply[Plus,matrow[n]]<10^(len-1),n++;
matrow[n]=Sum[matrow[k],{k,n-len,n-1}];];
While[First[matrow[n]]≤10^(len-1),res=list[res,matrow[n]];
n++;
matrow[n]=Sum[matrow[k],{k,n-len,n-1}];];
vecs=Apply[List,Flatten[res,Infinity,list]];
Map[(#-10^Range[len-1,0,-1])&,vecs]]

Code to solve a given repfigit equation. It gives out information as it runs.

SetSystemOptions["LatticeReduceOptions"->{"LatticeReduceRatioParameter"->.999}];

integerNullSpace[vec:{_Integer..}]:=Module[{mat,hnf},mat=Transpose[Join[{vec},IdentityMatrix[Length[vec]]]];
hnf=Last[HermiteDecomposition[mat]];
LatticeReduce[Map[Drop[#,1]&,Drop[hnf,1]]]]

SetOptions[LinearProgramming,Method->"CLP"];
keithSolutions[onulls_]:=Catch[Module[{nulls,len=Length[onulls],len2=Length[onulls[[1]]],order,stack,soln,solns={},badvar,val,counter=0,vval,eps=10^(-8),rndvar,bnds1,bnds2,obnds,mat,vec,b2b,gcds},
nulls=Reverse[onulls[[Ordering[Map[Norm,N[onulls]]]]]];
nulls=Transpose[nulls];
gcds=Apply[GCD,nulls,1];
nulls=nulls/gcds;
mat=N[Join[nulls,-nulls]];
bnds1=Join[ConstantArray[0,len2],Ceiling[ConstantArray[-9,len2]/gcds]];
bnds1[[1]]=1;
bnds2=Table[{-Infinity,Infinity},{len}];
While[bnds2=!=obnds,
obnds=bnds2;
Do[vec=UnitVector[len,j];
Quiet[soln=LinearProgramming[vec,mat,bnds1,bnds2]];
If[Head[soln]===LinearProgramming||!FreeQ[soln,Indeterminate],Throw[{0,{}}]];
val=vec.soln;
bnds2[[j,1]]=Ceiling[val-eps];
Quiet[soln=LinearProgramming[-vec,mat,bnds1,bnds2]];
If[Head[soln]===LinearProgramming||!FreeQ[soln,Indeterminate],Throw[{0,{}}]];
val=vec.soln;
bnds2[[j,2]]=Ceiling[val-eps];
,{j,len}];
];
order=Ordering[bnds2[[All,2]]-bnds2[[All,1]]];
bnds2=bnds2[[order]];
nulls=nulls[[All,order]];
mat=N[Join[nulls,-nulls]];
stack={{bnds1,bnds2},{}};
While[stack=!={},
counter++;
{bnds1,bnds2}=stack[[1]];
stack=stack[[2]];
rndvar=RandomInteger[{1,len}]*RandomChoice[{-1,1}];
vec=UnitVector[len,Abs[rndvar]];
Quiet[soln=LinearProgramming[vec*Sign[rndvar],mat,bnds1,bnds2]];
If[Head[soln]===LinearProgramming,Continue[]];
If[!FreeQ[soln,Indeterminate],Continue[]];
vval=Ceiling[soln.vec-eps];
bnds2[[Abs[rndvar],(1-Sign[rndvar])/2+1]]=vval;
badvar=Position[soln,(a_/;Chop[a-Round[a]]=!=0),{1},1,Heads->False];
If[badvar==={},
soln=Round[soln];
solns={Flatten[(gcds*nulls).soln],solns};
Print["sol ",Flatten[(gcds*nulls).soln]];
Do[
b2b=bnds2;
Do[b2b[[k]]={soln[[k]],soln[[k]]},{k,j-1}];
b2b[[j]]={bnds2[[j,1]],soln[[j]]-1};
stack={{bnds1,b2b},stack};
b2b[[j]]={soln[[j]]+1,bnds2[[j,2]]};
stack={{bnds1,b2b},stack};
,{j,Length[soln]}];
Continue[]
];
badvar=badvar[[1,1]];
val=soln[[badvar]];
b2b=bnds2;
b2b[[badvar]]={bnds2[[badvar,1]],Floor[val]};
stack={{bnds1,b2b},stack};
b2b[[badvar]]={Ceiling[val],bnds2[[badvar,2]]};
stack={{bnds1,b2b},stack};
];
Throw[{counter,Partition[Flatten[solns],len2]}]
]]

We generate all Keith number equations up to size 34.

Do[keqns[j]=keithEquations[j];Print[j];
Do[nulls[j,k]=integerNullSpace[keqns[j][[k]]],{k,Length[keqns[j]]}],{j,2,34}];

Here we show a solving run up to size 20. Takes around a minute.

Timing[knums=Table[Print["digits ",j];
Print[Timing[digits=Table[{j,k,First[Timing[knumlist=keithSolutions[nulls[j,k]];]],knumlist},{k,Length[keqns[j]]}]]];digits,{j,2,20}]]

Now find all repfigits in 21-24 digit range.

Timing[knums=Table[Print["digits ",j];
Print[Timing[digits=Table[{j,k,First[Timing[knumlist=keithSolutions[nulls[j,k]];]],knumlist},{k,Length[keqns[j]]}]]];digits,{j,21,24}]]

Now find all repfigits in 25-28 digit range. Takes maybe an hour (earlier versions took a couple of days at least).

Timing[knums=Table[Print["digits ",j];
Print[Timing[digits=Table[{j,k,First[Timing[knumlist=keithSolutions[nulls[j,k]];]],knumlist},{k,Length[keqns[j]]}]]];digits,{j,25,28}]]

Here I show the 33 and 34 digit runs. This is the largest done to date.

Timing[knums=Table[Print["digits ",j];
Print[Timing[digits=Table[{j,k,First[Timing[knumlist=keithSolutions[nulls[j,k]];]],knumlist},{k,Length[keqns[j]]}]]];digits,{j,33,34}]]